{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "大家好，这篇是有关台大机器学习课程作业三的详解，题目同Coursera。\n",
    "\n",
    "我的github地址：  \n",
    "https://github.com/Doraemonzzz\n",
    "\n",
    "个人主页：  \n",
    "http://doraemonzzz.com/\n",
    "\n",
    "作业地址:  \n",
    "https://www.csie.ntu.edu.tw/~htlin/course/ml15fall/\n",
    "\n",
    "参考资料:  \n",
    "https://blog.csdn.net/a1015553840/article/details/51085129  \n",
    "http://www.vynguyen.net/category/study/machine-learning/page/6/  \n",
    "http://book.caltech.edu/bookforum/index.php  \n",
    "http://beader.me/mlnotebook/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 1\n",
    "\n",
    "$$\n",
    "E_{\\mathcal D}[E_{in}(w_{lin})] = \\sigma^2(1-\\frac{d+1}N)\n",
    "$$\n",
    "\n",
    "这个结论的理论推导可以参考我写的Learning from data习题Exercise 3.3,3.4。\n",
    "\n",
    "对于此题来说，直接把$N$解出来即可\n",
    "$$\n",
    "N=\\frac {d+1}{1-\\frac{E_{\\mathcal D}[E_{in}(w_{lin})]}{\\sigma^2}}\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "44.99999999999996\n"
     ]
    }
   ],
   "source": [
    "def f(d,delta,Ein):\n",
    "    return (d+1)/(1-Ein/(delta**2))\n",
    "\n",
    "print(f(8,0.1,0.008))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "所以$N\\ge45$即可"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 2\n",
    "\n",
    "这题实际上是Learning from data Exercise 3.3以及Problem 3.10的结论，这里一并给出。\n",
    "\n",
    "(1)$H$是对称矩阵\n",
    "\n",
    "(2)$H^K=H$ ($K$为任意正整数)\n",
    "\n",
    "(3)$H$的特征值$\\in\\{0,1\\}$\n",
    "\n",
    "(4)$H$是半正定矩阵\n",
    "\n",
    "(5)$trace(H) = d +1$\n",
    "\n",
    "(6)$H$有$d+1$个特征值为$1$\n",
    "\n",
    "(1)\n",
    "$$\n",
    "\\begin{aligned}\n",
    "H^T&=(X(X^TX)^{-1}X^T)^T\n",
    "\\\\&=X((X^TX)^{-1})^TX^T\n",
    "\\\\&=X((X^TX)^T)^{-1}X^T\n",
    "\\\\&=X(X^TX)^{-1}X^T\n",
    "\\\\&=H\n",
    "\\end{aligned}\n",
    "$$\n",
    "(2)直接验证即可，先来看$K=2$的情形\n",
    "$$\n",
    "\\begin{aligned}\n",
    "H^2&=X(X^TX)^{-1}X^TX(X^TX)^{-1}X^T\n",
    "\\\\&=X(X^TX)^{-1}X^T\n",
    "\\\\&=H\n",
    "\\end{aligned}\n",
    "$$\n",
    "那么对于任意$K$\n",
    "$$\n",
    "\\begin{aligned}\n",
    "H^K&=H^2H^{K-2}\n",
    "\\\\&=HH^{K-2}\n",
    "\\\\&=H^{K-1}\n",
    "\\\\&=...\n",
    "\\\\&=H\n",
    "\\end{aligned}\n",
    "$$\n",
    "(3)因为$H^K=H$，所以对于​$H$的任意特征值$\\lambda$\n",
    "$$\n",
    "\\lambda^K=\\lambda 恒成立\n",
    "\\\\\\lambda=0或1\n",
    "$$\n",
    "(4)$H$为对称矩阵且特征值$\\in\\{0,1\\}$，所以由线性代数知识可知$H$半正定\n",
    "\n",
    "(5)利用迹(trace)的性质$trace(AB) = trace(BA)$\n",
    "$$\n",
    "\\begin{aligned}\n",
    "trace(H)&=trace(X(X^TX)^{-1}X^T)\n",
    "\\\\&=trace(X^TX(X^TX)^{-1})\n",
    "\\\\&=trace(I_{d+1})(注意H^TH为(d+1)\\times (d+1)阶矩阵)\n",
    "\\\\&=d+1\n",
    "\\end{aligned}\n",
    "$$\n",
    "(6)我们知道对称矩阵必然相似于对角阵，所以存在可逆矩阵$P$，使得$H=P^{-1}AP$，那么\n",
    "$$\n",
    "d+1=trace(H)=trace(P^{-1}AP)=tr(PP^{-1}A)=tr(A)\n",
    "$$\n",
    "而$A$为由$0,1$构成的对角阵(因为$A$和$H$相似且$H$的特征值$\\in\\{0,1\\}$)，所以$H$一共有$d+1$个特征值为$1$\n",
    "\n",
    "有了以上结论，看此题的选项就很轻松了。\n",
    "\n",
    "(a),(d),(e)成立，(c)错误，稍微要看一下的是(b)，我们知道$H$的特征值$\\in\\{0,1\\}$，所以$H$不一定可逆，因此(b)也错误。\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 3\n",
    "\n",
    "先对原式进行变形\n",
    "$$\n",
    "[\\![y\\neq sign(w^Tx)]\\!] \\Longleftrightarrow[\\![y^2\\neq ysign(w^Tx)]\\!]\\Longleftrightarrow[\\![1\\neq sign(yw^Tx)]\\!]\n",
    "$$\n",
    "令$s=yw^Tx$，所以几个误差分别可以写成\n",
    "$$\n",
    "err_1=[\\![sign(s) \\neq 1]\\!]\n",
    "\\\\err_2=max(0,1-s)\n",
    "\\\\err_3=(max(0,1-s))^2\n",
    "\\\\err_4=(max(0,-s))\n",
    "\\\\err_5=\\theta(-s)=\\frac{e^{-s}}{1+e^{-s}}=\\frac{1}{1+e^{s}}\n",
    "\\\\err_6=e^{-s}\n",
    "$$\n",
    "接着作图。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEGCAYAAABxfL6kAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzs3XdUVNfax/HvBoZeLDQVwa6o2Gts\nWMASNe1NuTdNb7ppxvRyU41Rb5Jriqbc9N6LGjXYiBorGqNGBCsKiiC9MzD7/eMMiAQQdWAAn89a\nLIczZ855xpjf7Nlnn72V1hohhBBNh4O9CxBCCGFbEuxCCNHESLALIUQTI8EuhBBNjAS7EEI0MRLs\nQgjRxEiwiwZDKRWjlOp6Hq8brZQadZZ9XJRSqpbHu0cp5aqUWqOUClNKPayU8lZKvaOUGlnLYxQr\npVxreH6mUiq00ra7lVIjanN8IWoiwS4akkKgCEApNVsptU8ptcH6c0ApNaNsR6VUf6XUNOuvw4Ah\n1R1UKdUM2A7EWD88Kv7kKaW0Uqp/hZc4AU8BJYAHcJ3WOhsYAyRWcfz/KaWmV9pcDJir2NdFKeVm\nPX6+UspZKdVCKdUGSAJuUkoFKqVCavpgEKImTvYuQAilVHuMYPYFLldKLcMIxYXAEutu12MEbZmb\ngR3Wx6VAfoXjOQNmbb37TmudCfSs8Hxr4EZgOLAU+Fhrfdj6nCtwHNBAcyAcWGdtqXtqrQ9Z93PR\nWhdZD1lViJdorUut+3YA+mutvwUigf8AFuBLjA+yL4CJ1m0AbwKOwNPA7rP+BQpRiQS7aAjKgrns\nTzOwARgBTKuw304ApZQf8H/AJUqpe4BAoFQpdbN1PxfgSmB/2QuVUl7AR8B9QAQQAlwNPAZkVTiH\nMzASGAQMAGKBdOAuwEEpFQMEAzlKqT5a6xyMD4G/UUp5WM93N/AGgNZ6iVJqKJADpGmt37Xu6w+k\nAO8BPwH/0VpLqIvzIl0xwu601ke11j8DGUAUkAzMAy6v9PO49SVvALO11gO01gMwAvv5st+11mFa\n64qh/jEwEFgEvIPR5ZMC/AgoYKRS6iPr7gXAIYxukfVAprWujsCb1vOtweieyanwNl5XSiUppbKU\nUjdat+3C6MrppbWeZ63FE5iK8U1kmnXbNOAq4Evrt4ynge8q98ELUVsS7MLulFKOSqmBGC3vZRgt\n5UCtdR/gn8B4YBzQQSnlgBHI71Y4RDPO7IpRSimXCs+bAIvWejVGK73MTK31M0Aup7+9BgGngEcw\nwv0Z4Kj1GGX98MHAwUpv4z6tdRuttY/W+lPrtl5a66e01ukV9vuH9RyfAD2UUtHAMYxWelfrN4vZ\nwAStdWxNf29CVEeCXTQENwFvAa7AlVrr9RWeew6YYH2sMbpQngI2l10ABW4F/lfh9xhgu1LKscLr\nyvq6l2O00v2ApUqpzpzZd98CeADjm8MlGC3rm6yPWymlmgNelcK6SlrrvCo2f4rxAXYlsFtrHW79\nwMnHaLVfAvhqrXee7fhCVEf62EVD8LnW+kOl1AaMbhIw+rrBuID5h/WxyXqRM6zshdbuil8wLkLe\nWk0gaox+8xsxuj+CgVSM/vVVwBPlO2q9XSk1HPgdGIzxbWGS1rpAKfUDRst6TXVvxPph0qyG96qs\n5x8ItFVK/Q/4DSPwfwd6Aa/X8Hohzkpa7MLutNbFVWy+XykVBHTBuKjqCsyquINSyhf4CrgHuAX4\nUSk1popjKYzRJeMxQv4GIENr/T1GuGdX2r81Rqt/BUaf/Hbr9h8wPiCWUAWl1FSMfvUpVTx3jVLK\nHbgTYxTMGOt55wFfaK1PYlwc7gN8XdXxhagtCXZhd9Y+9tYY3SCl1s1HgF8x+rr7ANFAL6WUyfqa\n8cBG4A2t9TKt9UbgOuADpdRSpdTUCv3sCiPID2EE+4/A2wBa6y8xRqhUdADjA8MFeBKjL/xKjMB9\nCHi30o1EXTBa2fdgfGv4CGOUTqC11o7AZ4C71vq/Wuv+wItAjtb6ANDbeoE3EGMI52ql1G3WvxMh\nzpl0xYiGQGN0cfwOJCulfgVaAvdqrdcAWMe2LwKuU0rNwhjzfaPWekv5QbTeopTqDjwIPI/xYVCE\ncRF0QaUbT++p8Lsn1u4ea6v6F4zW85Va62Sl1A3ADIyRMPFKqe3A00qpy6396F8Br2qtV1Q4/jvA\nGuvFXgvwb631qQrPOwPOSqkQjFE97wO3aK1LlFITrOeLxxhTL8Q5UbKCkmholFJelYYSVn6+GxCn\n5R+vEFWSYBdCiCZG+tiFEKKJkWAXQogmxi4XT319fXW7du3scWohhGi0tm/ffkpr7Xe2/ewS7O3a\ntSMmJsYepxZCiEZLKZVQm/2kK0YIIZoYCXYhhGhiJNiFEKKJkTtPhRBNmtlsJjExkcLCwrPv3EC4\nuroSFBSEyWQ6r9dLsAshmrTExES8vLxo164dtVzP3K601qSlpZGYmEj79u3P6xi16opRSgUopf6o\n4fn3lVKblFJPnVcVQghRRwoLC2nZsmWjCHUApRQtW7a8oG8Yte1jfxlwq6aIKwFHrfVQjBVuOp93\nNUIIUQcaS6iXudB6zxrs1vmt8zDWoaxKOPCN9XEUxsrvdSK3OJe5W+eSU1zt/FBCCHHRqzHYlVLO\nwL8xFiOojgfGtKhgrOYeUM2xbi9buiw1NfV8auVg1kG+3vc1T2x4Aou2nNcxhBDC3k6ePMmIESPO\nvuN5OluL/TFgkdY6s4Z9cjndTeNZ3TG11u+WrSLv53fWO2Kr1NuvNw8NfIjoY9G8t/u98zqGEELY\nU0ZGBjfffDN5eVUtiWsbZxsVMw4Yo5S6G+ijlHpPa31rpX22Y3S/bAZ6A3G2L/O0f3b7J7tSd/Hm\nH2/So2UPhrUZVpenE0I0Ic8t+Yu9xyuvhHhhurf25pkpPap9Pj8/n5tuuomUlBTCwsJ46aWX+Prr\nr7nssstsWkdFNbbYtdYjrauoh2OsKPOqUmp2pd1+Am5USr0KXIOx+kydUUrxzNBn6NS8E4+se4TE\nnMS6PJ0QQlyQd999l549e7Ju3TpOnDjBkSNH8PHxqdNz1nocuzXcAZ6qtD1bKRUORADztdZZNquu\nGu4mdxaEL+C6pdfxQPQDfDrxU1ydXOv6tEKIRq6mlnVdiYuLY+PGjURHR5OZmUlSUhK9evWq03Pa\nZEoBrXWG1vobrXV1I2dsLtg7mLkj57IvfR8vbH4BWQlKCNEQde3alZkzZxIdHc3s2bMJDg6u83M2\n6rliRgaN5K7ed7H44GK+ifvm7C8QQoh6dtttt7F8+XJGjhzJ22+/Tdu2bev8nI1+SoE7e9/JnlN7\nmLttLl1bdKWPfx97lySEEOU8PDz45pu/Nzyjo6Pr7JyNusUO4KAceGnESwS6B/Jg9IOcKjhl75KE\nEMKuGn2wA/i4+LBg9AKyi7N5+LeHKbGU2LskIYSwmyYR7ABdW3Tl6aFPE3Myhv9u/6+9yxFCCLtp\nMsEOMKXjFP7R7R98svcTlh9ebu9yhBDCLppUsAM8POBh+vn34+nfnyY2Ldbe5QghRL1rcsFucjTx\nSvgr+Lj4cP/a+0krSLN3SUIIUS4rK4uJEycSGRnJFVdcQXFxsc3P0eSCHcDXzZfXRr9GemE6D/72\nIGaL2d4lCSEEAJ9//jmzZs0iKiqKwMBAVqxYYfNzNPpx7NXp4duDZy95lsfXP878rfN5csiT9i5J\nCGFvyx+D5N22PWZgGEycW+3TlScBW7hwYflzqamp+Pv727YemmiLvczkDpOZ1mMaX8V9xffx39u7\nHCHERajyJGC7du0CYNOmTWRkZDBkyBCbn7PJttjLzOw3k/iMeGZvmU3HZh3lzlQhLmY1tKzrSlWT\ngAUFBXHvvffy/fd10+Bs0i12AEcHR+aPnE8rj1bMXDuT5Lx6m6dMCCH+NglY69atufrqq3nppZcI\nCQmpk3M2+WAH487U10e/TkFJAQ+sfYCi0iJ7lySEuEhUngQsOjqaHTt28OKLLxIeHs7XX39t83Mq\ne0x3O2DAAB0TE1Pv511zdA33r72fKR2m8OLwFxvdyuVCiHMXGxtLaGiovcs4Z1XVrZTarrUecLbX\n1qrFrpRqoZSKUEr5nmeNDcKY4DHM6DODJYeW8OneT+1djhBC1ImzBrtSqjmwFBgErFVK/W0laqWU\nk1LqqFIq2voTVge12sQdve5gbPBYXtn+ChuPb7R3OUIIYXO1abH3AmZprV8EfgX6VbPPl2Xro2qt\nbTxQ1HYclAMvDn+RDj4deOi3hziSdcTeJQkhhE2dNdi11r9prTcrpUZitNo3VbHbEGCyUmqrUup9\npVSDHkbpYfLgjTFv4KScuGfNPWQV1fkyrUIIUW9q28eugGuBDKCq+/O3AeO01oMAEzCpimPcrpSK\nUUrFpKamXkDJthHkFcSC0Qs4nnucB6Nl2gEhRP1JT09n5cqVnDpVNwsD1SrYteFuYBcwtYpddmmt\nT1gfxwCdqzjGu1rrAVrrAX5+f+umt4t+Af14ZugzbEnewktbXpIFsYUQdS4jI4PJkyezdetWRo8e\nTV00dGtz8fRRpdRN1l+bAZlV7PapUqq3UsoRuBz404Y11qnLOl3Gv3r+i2/jv+WLfV/YuxwhRBO3\na9cuXn31VZ588knGjx/Pjh07bH6O2vSFvwt8o5S6FdgDJCqlZmutn6qwz/PAF4ACFmutV9m80jp0\nf7/7OZx1mPnb5hPiHcLwNsPtXZIQog7M2zqPfen7bHrMbi268eigR6t9vrpJwNatW8fWrVt5+umn\nbVoP1O7iaYbWOkJrPVJrPUNr/VelUEdrvUdr3UtrHaa1bnTTKDooB+aOmEvnZp15+LeHOZR5yN4l\nCSGaiKomAdNa8/XXX9O8eXNMJpPNz9mgR6/UJ3eTO2+MeYN//PIP7l59N19e+iXNXJvZuywhhA3V\n1LKuK1VNAtarVy8WLlzIv//9bxYvXsy1115r03NeFHPF1FYrz1a8NuY1UvJTeCD6AcylMlJGCHFh\nKk8CtnXrVj755BMAMjMzadbM9g1ICfZKevv15vlhzxNzMobZW2bLSBkhxAWpPAnYtGnT+PTTTxk5\nciSlpaVERkba/JzSFVOFSztcyqGsQ7y76106+HTg5h4327skIUQj5eHhwTfffHPGtpUrV9bpOSXY\nq3F3n7s5nHWYV2JeIcgriLHBY+1dkhBC1Ip0xVSjbE6ZMN8wHlv3GLtTG+z0N0IIcQYJ9hq4Obnx\n+pjXaenWknvW3ENSbpK9SxJCiLOSYD+Llm4tWTRuEWaLmRmrZpBdnG3vkoQQokYS7LXQwacDr41+\njaM5R3lgrQyDFEI0bBLstTQwcCDPX/I8W5O38uymZ2UYpBDigpw8eZK+ffvWybFlVMw5mNJxCom5\niSzauYi2Xm25s/ed9i5JCNFIPfTQQxQUFNTJsSXYz9Gdve4kMSeRhTsX0sazDVM6TrF3SUKIWkqe\nM4eiWNtOAuYS2o3AJ56o9vmqJgFbs2YNHh4eBAYG2rSWMtIVc46UUjw79FkGBQ7i6Y1Psy15m71L\nEkI0YFVNAvbCCy8wd+7cOjuntNjPg8nRxKvhr3LT8puYuXYmn036jPY+7e1dlhDiLGpqWdeVypOA\n/fjjj8yYMaNO5ogpIy328+Tj4sOicYtwcnDirlV3caqgbpa4EkI0bpUnAXvzzTdZuHAh4eHh7Ny5\nk1tvvdXm55RgvwBtPNuwcOxC0gvTmbFqBnnmPHuXJIRoYCpPAnbw4EGio6OJjo6mT58+vPfeezY/\nZ20Xs26hlIpQSvnavIJGrqdvT14e9TLxGfHMip4lY9yFEGcomwRs3bp1LF26FG9v7/LnoqOj6+Sc\ntVnztDmwFBgErFVKVbkStVLqfaXUJqXUU1U935SNDBrJM0OfYePxjTyz8RkZ4y6EsKvaXDztBczS\nWm+2hnw/4NeKOyilrgQctdZDlVIfKKU6a63310G9DdYVna8gJT+FN3e+iZ+7Hw/0f8DeJQkhLlJn\nDXat9W8ASqmRGK3256vYLRwom3A4ChgOXFTBDnB7r9tJyU/hgz0f4O/uz/Wh19u7JCHERahWwx2V\nUgq4FsgAqupE9gDKpj5Mx2jVVz7G7cDtAMHBwedTa4OnlOKJwU9wquAU87bOw9fNl/Htxtu7LCHE\nRaZWF0+14W5gFzC1il1yATfrY8+qjqu1fldrPUBrPcDPr8pu+ibB0cGReSPn0ce/D4+vf1xuYBJC\n1LvaXDx9VCl1k/XXZkBmFbttx+h+AegNHLFJdY2Uq5Mrb4x5g7Zebbl/zf3EZ8TbuyQhRANRUlJC\ncHAw4eHhhIeHs3u37RfxqU2L/V3gRqXUOsARSFRKza60z0/WfV4FrgF+sW2ZjY+Piw9vj3sbNyc3\n7lp1F8l5yfYuSQjRAOzatYt//OMf5WPZw8LCbH6O2lw8zQAiKm1+qtI+2UqpcOt+87XWWTarsBFr\n5dmKtyLe4ublN3Pnyjv5eOLH+Lj42LssIS5a67+J59SxXJse07etJyOu6VLt85UnAevRowdLly5l\n7dq1hIWF8c477+DkZNvZXWx256nWOkNr/Y3WWpqmFXRp3oXXx7zO0ZyjzFg9g3xzvr1LEkLUo8qT\ngPXo0YNVq1axdetWzGYzy5Yts/k5ZRKwejAwcCD/GfUfZkXP4oHoB3hjzBs4OzrbuywhLjo1tazr\nSuVJwDIzM2nVqhUAAwYMYP9+248Ml7li6snY4LE8O/RZNh7fyOPrH6fUUmrvkoQQ9aDyJGCfffYZ\nf/75J6Wlpfz000/07t3b5ueUFns9uqLzFWQXZ/NyzMt4b/Hm6SFPY9wiIIRoqm677TamT5/Ohx9+\niLe3N3PmzOGGG25Aa83UqVMZN26czc8pwV7Pbu5xM1lFWfxv9//wcfZhZv+Z9i5JCFGHyiYBq2jX\nrl11ek4Jdju4t++9ZBZl8v6e92nm0oxpPafZuyQhRBMiwW4HSimeHPwkOcU5vLL9FXxcfLii8xX2\nLksI0UQ0uounSXGxTWJaXEcHR+YMn8OwNsN4dtOzrEpYZe+ShBBNRKMK9oRdO/nq6YdZ/f5bWEob\n/6gSk6OJV0e9Si/fXjyy7hE2n9hs75KEEE1Aowr24J69GDj1Kv5cuYwf5z9PUX7jv9nH3eTOm2Pf\npJ1PO+5bcx+7Uuv2oooQoulrVMGuHBwYef10Im6/h4Rdf/DVM4+QfSrF3mVdMB8XH94Z9w6+br7c\nuepOYtNi7V2SEKKOzZgxgyVLltTJsRtVsJfpNXYCVz7+HNmpKXzx5IMkH2z8a3r4ufvxXuR7eJo8\nuWPlHRzIOGDvkoQQdWT9+vUkJyczZcqUOjm+sseFyAEDBuiYmJgLPs6pYwn8OO858rOymHTfQ3Qe\nONQG1dnX0eyjTFsxDY3mowkfEeIdYu+ShGjUYmNjCQ0NBWDtR++SknDIpsf3D+nA6Gm3V/t85UnA\nFixYQFhYGJMmTWLUqFFcdtllZ627jFJqu9Z6wNlqapQt9jK+bUP45+xX8A0OYfErc4hZ8kOjHzET\n7B3Me5HvYdEWbo26laTcpLO/SAjRYFWeBGzRokV0796dRx55hK1bt/LGG2/Y/JyNusVexlxUyPKF\nr7J/y0Z6jZvAmOl34mjjaTDrW1x6HNN/nY6Psw8fTfiIAI8Ae5ckRKNUVcu3Pt11111s3LiR5s2b\nk5mZSU5ODgsXLmTChAnExsby5JNP8sMPP/ztdRdti72MycWVKTMfY+Bl/8euVSv47sWnyM9u3FPC\nd23RlXfGvUNGUQa3Rt3KqYJT9i5JCHEeKk8C9q9//YtDh4zuoJiYGEJCbN/detYWu1LKB/gKY/Wk\nPOBarXVxpX2cgEPWH4B7tdbVrvdk6xZ7RXvXryXqndfxaNacyx56Cv92HerkPPVlx8kd3LnqToK8\ngvhw/IeyUIcQ58jeLfa8vDymT59OcnIy3t7efP7559x6662cPHkSs9nMd999R5s2bf72ugtpsdcm\n2GcA+7XWK5VSbwHLtdaLK+3TDyPwHz3726zbYAdIPhDPz6+8SGFeLhPvnkWXwcPq7Fz1YdPxTdyz\n+h46Ne/Ee5Hv4eXsZe+ShGg07B3s56tOu2K01ou01iutv/oBVQ0cHwJMVkptVUq9b23B201gpy5c\nP+e/+IW0Z8mrL/H7N5+hLRZ7lnRBhrYeyn9H/5f4jHhmrJJVmIQQNat1H7tSaijQXGtd1X3v24Bx\nWutBgAmYVMXrb1dKxSilYlJTU8+74NrybN6Ca55+iZ6jI9j8/VcsfnUOxQWNNxBHBo1k/sj57Dq1\nS5bYE+IcNbbRchdab62CXSnVAngD+Fc1u+zSWp+wPo4BOlfeQWv9rtZ6gNZ6gJ+f33kVe66cTCYi\n77iP0dNu5+D2rXzx1ENkJp84+wsbqIiQCOaOmMsfKX9w9+q7JdyFqAVXV1fS0tIaTbhrrUlLS8PV\n1fW8j1GbPnZnYDkwt0KXTOV9vgFeBPYAK4E5Wutqpyus6z72qiTs3snSBfPQFgsT75lFx/6D6/X8\ntrTs0DIe3/A4/fz7sXDsQtxN7vYuSYgGy2w2k5iYSGFhob1LqTVXV1eCgoIwmUxnbLflxdO7gDnA\nn9ZNawGT1vqpCvv0BL4AFLBYa/1kTce0R7ADZKUks/jVl0g5fJDBV1zDJddcj4ODY73XYQtl4d4/\noD9vjnlTwl2Ii4DNgr0u2CvYAUqKi1nz4dvsXhNFcM/eXHr/I7h7N84hhL8c+oUnNjwh4S7EReKi\nukHpXDg5OxN5x31E3nkfSXF7+fSx+zkev8/eZZ2XSztcypzhc9h+cjv3rrlX+tyFEMBFGOxlwkZH\n8o8XXsbR0ZGvn32MP35d2mgurlR0aYdLeXH4i8ScjOHeNfdSUFJg75KEEHZ20QY7QED7jtzw0mu0\n692XNR+8zbI3Xqa4sPEF4+QOk8vD/Z7V90i4C3GRu6iDHcDV05PLH/43w669kX0b1/HZ4w+QcsS2\n03rWh8kdJjN72GwJdyGEBDsYKzMNufJarn7qRYoL8vniqQcbZdfMlI5TysNdxrkLcfGSYK8guGcv\nbpr3OsE9erHmg7dZ/MocCnNz7V3WOZnScQpzhs9hx8kd3LHyDnKKc+xdkhCinkmwV+Lu04wrHn2G\nkTf8i0M7tvLJo/dyPL5xrUF6aYdL+c+o/7AnbQ+3Rt1KZmGmvUsSQtQjCfYqKAcHBk65kuuen4+D\ngwNfPfMoW376tlFNJBYREsFro1/jQMYB/hX1L9IK0uxdkhCinkiw16BVp67cOO91Og+6hA1ffsz3\nLz1DbnrjCciRQSNZOG4hiTmJTP91OifzTtq7JCFEPZBgPwsXdw8mz3yUiNvuIWnfXj5++B7it/xu\n77JqbUirIbw17i1S8lOYtmIax3OP27skIUQdk2CvBaUUvcZN4MZ5r+HjH8iSV19ixaIFFOU3jlEn\n/QP6827Eu2QVZ3Hzips5mn3U3iUJIeqQBPs5aNE6iH+88B+GXHkte9et4ZNH7iVx31/2LqtWevn1\n4oPxH1BUUsS0FdM4lNn4xuoLIWpHgv0cOTo5MezaG7n2uXkoBd88+zgbvvqE0hKzvUs7q24tuvHB\n+A/QaKb/Op249Dh7lySEqAMS7OepTddQbpr/Bt1HjWHLj9/w5b8fJi3pmL3LOqtOzTvx0YSPcHZ0\nZvqv09mZstPeJQkhbEyC/QI4u7kz4a6ZTJ31BFmpKXz66H1sW/w9FkupvUurUYh3CB9P+JgWri24\nLeo2NiRtsHdJQggbkmC3gc6DL2Haywtp36c/6z7/kK+efqTBt95be7bm4wkf096nPfeuvpflh5fb\nuyQhhI1IsNuIR7PmTH3wSSbd9zAZJ443itZ7S7eWvD/+fXr79+bRdY/y9b6v7V2SEMIGzhrsSikf\npdRypVSUUupH6xqoVe33vlJqk1LqqaqevxgopQgdNoppryxqNK13L2cv3h73NqPajmL2ltm8/efb\njW7yMyHEmWrTYr8eeFVrHQkkAxMq76CUuhJw1FoPBToopTrbtszGpcrW+5IfsJQ2zNa7q5Mr/w3/\nL1M7TmXhzoXM3zYfi2480ycIIc7kdLYdtNaLKvzqB6RUsVs48I31cRQwHNh/ocX9TV4a/HAbhE6G\nbpPB09/mp7CVstZ7cI9erHpvIes++4C4jeuIvOM+/Nt1sHd5f+Pk4MQLw17A29mbz2I/I6soi+eG\nPYfJwXT2FwshGpRa97ErpYYCzbXWm6t42gNIsj5OBwKqeP3tSqkYpVRMamrqeRVLZgJkHIalD8DL\nXeDDSbD5bchKPL/j1YOy1vvkmY+Rk3aKzx6fybrPP8RcVGjv0v7GQTnwyMBHuLfvvSw5tIRZa2dR\nWNLw6hRC1EzVpj9VKdUCoyV+ldY6oYrnXwO+1FpvtnbLdNNaz6nueAMGDNAxMTHnV7HWkLIX9i6G\n2MXGY4A2/SF0KnSfCi0aXosYoDA3l3VffMju1b/i4x/AuNvuoV2vvvYuq0rfxH3D7M2z6evfl9fH\nvI6Pi4+9SxLioqeU2q61HnDW/c4W7NaLpcuBuVrrldXscxPgr7V+WSn1HBCntf6iumNeULBXduqA\nEfCxi+H4H8a2gJ6nQ96vGyhlm3PZyLG9u1n57ptknEii+4jRjLrpVty9G15wrji8gsc3PE4773a8\nNe4tAj0C7V2SEBc1Wwb7XcAc4E/rprWASWv9VIV9vIH1wGpgIjBEa51V3TFtGuwVZR6F2KVGyB/d\nDGho2el0yLfq02BCvqS4mC0/fcPWn77D2d2d8BtvofvIMagGUl+ZrSe2cv/a+3F3cmfRuEV0bdHV\n3iUJcdGyWbCfwwmbAxHAOq11ck371lmwV5STDPuWQuwSOLwedCn4BEPoFCPkgwaBg/2H8Z86lsDK\nd9/keHwsbbp1Z+wtM/ALbmfvss4QnxHPXavuIs+cx4LRCxjSaoi9SxLiolTvwX4u6iXYK8pPh7hl\nRsgfXAOlxeAZYIys6T4VQoaD41kHCNUZbbGwJ3oV6774iKK8XPpNnMLQ/7seF3d3u9VUWXJeMnet\nuosj2Ud4YdgLTO4w2d4lCXHRkWCvTmE27I+CvT/DgVVgzge3FtBtktFl0yEcnFzsUlpBbg4bvvyY\nXat/xaNZc0bdeAvdLhnZYLoGFw5YAAAgAElEQVRnsouzmbl2JtuSt3F/v/u5pectDaY2IS4GEuy1\nUZwPB1cbI2ziV0BRNrh4Q5fxRpdNp3Hg7FHvZSUfiGfV+29x8tB+2vboxdh/3UnLoOB6r6MqxaXF\nPLXhKZYfWc61Xa/l8UGP4+jgaO+yhLgoSLCfq5JiOPyb0ZLf9wsUpIOTG3QeB6GXQZdIcK2/kSsW\nSym7V//Khi8/obiwgL4TpjDkqutw9fCstxqqrU1bWLBjAR/u+ZDRbUczb+Q83Jzc7F2WEE2eBPuF\nKC2BoxutY+WXQG4yODob3TShU6HrJPBoWS+l5Gdnsf6Lj9kTvRI3Ty+GXXsjYWMjcWgAreQvYr9g\n7ta5hPmF8caYN2jh2sLeJQnRpEmw24rFAkkxRks+drExpFI5QrthRsiHTgGvuh/fffLQAdZ+/D+S\n9v2FX3A7wm++neCever8vGezKmEVj61/DF83XxaNXUSHZg3z5jAhmgIJ9rqgNZz402jFxy6GU/GA\ngraDTod885A6PL0mfvPvrPv8A7JTU+g0cCijbryFZgH2vXFoV+ou7ltzH8WlxbwS/gpDWw+1az1C\nNFUS7PUhZZ815H+G5N3Gtla9rTdEXQa+dTPJpbm4iO1Lf2LrT99iKS2h36WXM/jyq3Fxr/8LvWWO\n5x7n7tV3czjrME8OeZKru1xtt1qEaKok2Otb+iFryC+BxG3GNr9up+96Dehp87tec9PTWP/lx+xd\ntwZXL2+GXHEtvSMn4WSyz4yMucW5PLzuYTYkbeDm7jfzQP8HZMSMEDYkwW5PWUnGXa97FxsXYbUF\nmre33vV6mTFhmQ1D/uShA6z7/EOO7vkTH/8Ahl17ozH+3Q531pZYSpi/bT5f7vuS8LbhzBsxD3dT\nw7nRSojGTIK9ochNhbhfjJA//BtYSsC7jRHyoVMgeCjYoFWrtSbhzx2s++IjUhMO49++IyOvn05I\nWB8bvIlz93ns58zfNp8uzbvwxpg3ZAIxIWxAgr0hKsg0boTau9i4MaqkEDz8oNulRsi3HwWOF9aN\noi0WYjdEs+HrT8k5lUq73v0Y8c9pdlncY13iOh7+7WE8TB68Nvo1wvzC6r0GIZoSCfaGrigXDqw0\nQn5/FBTnGjdAdbVObdBxNJjO/6afkuJidv66lC0/fkNhfh7dh4dzyTXX4+Nfvy3n+Ix47ltzH6n5\nqTx7ybNM6TilXs8vRFMiwd6YmAvh0Foj5OOWQWEmmDyMu11Dp0LnCHDxOq9DF+bmsvXnb9mxfDHa\nYqHn6AgGX3Et3r5+Nn4T1csozODB3x5kW/I2pvWYxsx+M+WiqhDnQYK9sSo1w5H1RsjvWwp5qeDo\nAp3GWu96nQBuzc/5sDnpp9jy47fsXv0rSkGvcRMZdPnVeDavn7tFzRYz87bO4+u4rxnWZhjzR87H\n29m7Xs4tRFMhwd4UWEqNBUPKhlFmJ4KDE7QfaYR8t8ngeW4t7+zUFDb/8BV7olfh6GSid+QkBl32\nf/W2gtO38d8yZ/McgryCeH3M67T3aV8v5xWiKbBpsCulAoDvtNYjqnm+DbAFOGDddLXWutoVqyXY\nz4PWkLTDuBlq72JjUW/lYIyqKbvr1adNrQ+XmXyCTd9/Sez6aJycnek7cQoDplyJm+f5dfmci5jk\nGGZFz6LEUsK8kfMYEVTlPyshRCW2XBqvOfAlxpqm/arZ50ogQGv9Vm2Kk2C/QFrDyb+sa70uuaAF\nvdOSjrHpuy+J27QeZ1c3+k2cQt+JU+u8BX889zj3rbmP/Zn7ub/f/UzvMV3mdhfiLGwZ7N6AAn7W\nWodXs898YKx1vxVa6ydqOqYEu42dOnC6JX9ip7EtIOz0MoC1WNA79egRNn33Bfu3bMTk4kqviIkM\nmHxFnfbB55vzeXrj0/x65FfGBY/jhWEv4Ols/2mJhWiobN7HrpSKriHYRwMxQD6wCrhfa72rumNJ\nsNehzKNGK37vYji2BWNB786nQ/4sC3qnJR5ly0/fsm/Dbzg4OdIzPIKBU6/Cxz+gTsrVWvPJ3k/4\n7/b/EuwdzILwBTJDpBDVqO9gd9FaF1kfvwr8rrX+vtI+twO3AwQHB/dPSEio1XnFBShb0HvvYjiy\nwVjQu1nw6T75Ghb0zkw+wdbF3/FX9GpAEzp8NIMuv5oWrWvfj38utp7YysPrHqawpJDZw2cTERJR\nJ+cRojGr72CPBv4BZAFbgau01nHVHUta7HZQtqD33sXGmPnSYvAMhNDJRtCHDKtyQe/sU6nELP2B\n3at+paTETNchwxl0+dV1cidrcl4ys6JnsfvUbqb3mM59/e7DycF+i4wL0dDUWbArpcYA3bXWb1Z4\nbjTwFlAMvFvxuapIsNtZjQt6XwYdRv1tQe+8zAy2L/uZP6N+obiggJBefRkw+QpCevW16UXP4tJi\n5m6dy7fx3zI4cDDzR82XlZmEsJJx7KJ2alzQe6p1Qe/TszMW5uby56rl/LFiCXkZ6fgGt2PA5Cvo\nNmwkjk62my74x/0/MnvzbJq7NuflUS/Tx98+k5kJ0ZBIsItzV1IEh34zhlFWuaD3eHA17hYtMZvZ\n9/tvxCz5gbTEo3i2aEm/iVPpNW6CzRb8+CvtLx6MfpCTeSeZ2X8mN3W/SYZEiouaBLu4MKUlkPC7\ndaz80r8v6N3tUnBvgdaaIzu3E7P0B47u2YWzmxthYyfQb+IUvH39L7iM7OJsnv79aVYfXU1423Bm\nD5uNj0v93CUrREMjwS5sx2IxVoWKXWx02WRVvaD3yUMHiFn6I3Gb1gPQaeAQ+k2YSpvQHhfU0tZa\n83ns57yy/RX83fx5edTLMgWwuChJsIu6Ub6gtzXk0/ZjLOg9uHzxkOwSN3ZG/cLu1b9SmJeLX3A7\n+k6cSrfhozA5u5z1FNXZnbqbh357iJSCFGb1n8UNoTdI14y4qEiwi/qRsu90yJ8sW9C7D4ROwdxp\nIrGxifyxYgmnjh7B1dOLsLHj6RM56by7abKKsnjq96eIPhbN2OCxPD/seZklUlw0JNhF/Stb0Hvv\nYkiy/vf1C0V3m0yiSxh/bNnDgW1bAKObps/4ybTtEXbOre6yu1UXbF9AgEcAc0fMlVEz4qIgwS7s\nKyvRGFlTaUHv7Lbj2ZnSjN3bdlGYm0Pz1kH0GjueHqPG4uZ1bi3vP1P/5NF1j5Kcl8ydve/ktrDb\nZAEP0aRJsIuGo4oFvc2eQcS7DOfPo3DiaBKOJhNdBg+j17gJtOlW+4utOcU5zN48m2WHl9HPvx9z\nR8yllWerOn5DQtiHBLtomAoyIG6F0WVzYBWUFpGqgvizpDexCYUUFxXTMiiYXmPH033kWFw9azfb\n45KDS5i9eTaODo48O/RZIttF1vEbEaL+SbCLhq8o15jaIHYJ7I/CXJjPvoJgduV2IDndjJPJmS5D\nhxM2OrJWQyaPZR/j0fWPsvvUbq7qfBWPDHwEd5N7ja8RojGRYBeNi7kQDq4xRtjELeNkZgm7s9sS\nm+VPcYmmmX8APcIj6D5qbI0LcZstZhb+sZAP9nxAiHcI80bOo3vL7vX4RoSoOxLsovEqNcPhdRC7\nGPNfy4hPhr+yW3Esz7i4GtKjJz3GTqTTwCHVjovfcmILT6x/gvTCdO7sfSe3hN0iM0WKRk+CXTQN\n5Qt6Lybzj+X8laTZmxVAttkVFxcT3YYOp0fEZAI7dvlbV01WURZztsxh2eFlhPmG8eLwF2XxbNGo\nSbCLpsdigeN/oPf+xNFNUfyVqNmf05IS7Ujz5h6EjhhL6NgpNAs8c1TMiiMreGHTCxSXFvNA/we4\nrtt1OKiqFxgRoiGTYBdNm3VB76KdPxC3YTWxSaUk5jcDoFWgD6GjIug67oryRblT8lN4ZuMzbEja\nwJBWQ3hh2AsEegTa8x0Icc4k2MXF5dR+srd+w74N0exLNJNa5IlC0y64BaHhE+g05kqcXF35Nv5b\nXo55GSflxOODH2dyh8ky34xoNCTYxcUrI4FTG74k9vffiD1mJqfEFScHC506+BE6ehIOfQfw7y3P\nsjN1J+Ftw/n3kH/j737hUwwLUdck2IUAdNZxklZ/QuzGDcQnmSm0mHBxLKVDe1+SuvnxP70OJydn\nHh74MJd3ulxa76JBs2mwK6UCgO+01iOqed4E/AC0AN7XWn9Q0/Ek2IU9lGad5MivH7F/60YOJBVR\nZHHC5FhKSqtStgdkEhLWm6dHPEdrz9b2LlWIKtks2JVSzYEvAX+tdb9q9pkFeGutn1VKLQOu1Vrn\nVHfM8w327LxC4hZHUdJnAFzAvN5CqPx0SrZ9S9ruLRxLNkK+2MnCcf8C+vfpzT//7zmcPWSlJtGw\n2DLYvQEF/Ky1Dq9mn8XAY1rrvUqpx4AtWuu1lfa5HbgdIDg4uH9CQkKt3khFu5auwfTQ3eQ7ubAt\noBu/tw5jW0AohU4S8uL8uet8xhf9ToeCOApywbHEkVJHC+1au9Br0Ejaj78RFx9fe5cphO372JVS\n0TUE+2rgSq11ljXAs7XWX1V3rPNusefkE7fiN0wbf8O0aR0OmRlokzMl/QZRPGwUJYOGoc9x6lch\nyjy35C98XUqY6P0Lm3bsoFWKK67FjjgoC239nOnUtz8dJ9yAV2u5yUnYR22D3Vb3WOcCbkAW4Gn9\n3ea8vdwZePVEuHoiurSUgh07yF65kpyVqzBt2QBOTngMGYJXRARe48bi1LJlXZQhmihfTxfyix25\n/vpXiLgihXmbXmLnrvWEpXrjklRCwq9bWP3rFgJ8HOgY1p1Okdfh26W3XHAVDY6tWuxPA3u11t8p\npT4G3tFab6zuWLa+eKq1pnD3bnKiosiOWon56FFwcMC9f38j5CMjMAXKzSiiZrd8tI3k7EJ+ue/0\nGIF1iet4cfOLHM89zv+5D2b04VIS4xM5kesMKHzcoGNoBzqMvoI2fYfjZDLZ7w2IJq/OumKUUmOA\n7lrrNys8FwIsA1YBlwBDtNal1R2rLkfFaK0pio8n59coclZGUbT/AACuvXvhHRGBV2QkzsHBdXJu\n0bjN+Hw7+0/msnLWqDO255vzeXvX23zy1yf4uPjw0ICHGF3szaHV33Jwzz4SMp0p1Q6YHDXBIQF0\nGDqOdpdE1DgLpRDno97HsSulWgPDgV+11lk17Vufwx2LDh0mZ+VKcqKiKPzrLwBcunXDKzIC74gI\nnDt1kq/SAoAHvt7J9oQM1j0yusrn49LjeH7z8+xK3cXgVoN5cvCTtPduhzlxF0dXf87hP//gUIom\np8QVAN+WHrTvN5j2l0TQuksojk4yu6S4MHKDUhXMSUlGn3zUSgr++AO0xrl9e7wiI/GKjMC1e3cJ\n+YvYY9/vYm1cClueGFftPhZt4bv471iwfQEFpQXcGHojd/S+Aw+TBwD61AHSN37Foa3rOZyUR1K+\nNxYccHZ2pF337rQbPJp2vfvh1VJG2YhzJ8F+FuaUFHJXryY7Kor8rdugtBRTmzbWPvlI3Pr0RjnI\nDIAXk2d+3sPPfx5n59NnX1YvrSCN13a8xo8HfsTPzY9ZA2ZxaftLz2wYZCVStPNHjm5azuHDKRzO\nbU5uiTE0t4W/L8F9hxDSqy9tu4fh4i4rPYmzk2A/ByUZGeSuWUtOVBR5GzeizWac/PzwihiHV2Qk\n7gMGoORrdJM3Z1ksn25KIPaFCbV+za7UXby05SX2pO2hn38/Hh/8ON1adPv7jrkp6NilnNq6mIS4\nAyTkepNY0IwSiwPKQdGqczdCwvoQEtaXwE5dpNtGVEmC/TyV5uSQG/0bOStXkrt+PbqgAMdmzfAc\nOwbvyEjchw7FwdnZ3mWKOvBKVBwL1x7g4JxJ59QlZ9EWfjrwEwu2LyCrOIuru1zNvX3vxcelmjtX\nrQt6l+z5meO7t3E0x4OEAj+S890AcHZzo22PXgT37E3bHr3wDQqWb48CkGC3CUtBAbnr15OzchW5\na9diyc3FwdMTz9Gj8YoYh+eIETi4udm7TGEjb67Zz8tR8ex/cSImx3MP0qyiLBbtXMRXcV/h5ezF\nXb3v4pqu12ByqGEIZPmC3osp2LuaY1nOJBT6k1DoT1aeMbDM1dOLoNCetO0RRtvuYfi2DZGgv0hJ\nsNuYpbiY/E2byI6KInf1GkozM1GurniOGIFXZCSeo8Nx9PS0d5niAry77iBzlu1jz3Pj8XQ5/66Q\nuPQ4/hPzH7ac2EI773Y80P8BRrcdffZvAeYCOLi2fEHvrOxCEosDOObYlcRMJ7Iyjfv+jKDvQdvu\nYQR1D8MvuJ0E/UWiSQa71rpBjFrRJSXkx8SQExVFzspVlKSmokwmPC65BK/ICDzHjMGpeXN7lynO\n0ccbj/DM4r/Y/tQ4Wnpe2PxDWmvWJ63n5ZiXOZx1mAEBA3ho4EP0aNmjdgcoKYYj6yB2CcQuhfxT\nZFu8OeYxmGMlrUk8kUtWagoArh6eBHXvSVBoGG26huLXroP00TdRTTLYU4/m8P387bi4O+HiYcLV\n3cl47G6y/ml97HHmNlfrNieTo83fi7ZYKNj5pxHyUVGYjx8HR0fcBw3EOzISz7FjMfnLIg6NwVdb\nj/LYD7vZ+NgYWjezTReb2WLm+/jvWbRzERlFGUzpMIX7+t13bsvyWUrh6CZryC+B7CRwcCI7cCSJ\nrn04luVK4v4DZJ48AYCTiwutOnahddfutOkaSqsu3XD1kG+TTUGTDPbstAL2RCdRlG+mKL+EwvyS\n8sdF+SUUF5TU+HpHk0N5+Nf0oeBavq3sQ6J2Hwpaawr37iUnyrghqvjwYVAKt759y2+IMrVpc87v\nW9SPH/9I5IGv/2TtQ+G09/Ww6bFzinN4b/d7fLb3M5RS3Nj9Rqb3nI638zlOWmexwPEdsPdno8sm\n4wgoBwi+hJy2ERxX7Tl+7CRJcbGkHDmItlhAKXyDgmndNZQ2XbvTuksoPgGBDeLbrzg3TTLYz8Zi\n0RTnl1BYHvanQ78o30xRXkmFDwXrc9ZtxYXVzoAAVPOh4GEq/zbg6mHC1dP43dXT+GBQKccoWLuK\n3JUrKdq3DwDXHj3Kb4hyaS+zBDYky3afYMbnO1gxcwTdAutmltCk3CRe2/Eayw8vx9vZm3/1/Bf/\nDP0nbk7n8Q1Bazi5B/YuNkI+1fg3RpsB0H0q5g7jOZFWzPG4vSTFx3Iifh9F+XkAeDRrTusuobTq\n3JXATl0I6NAJZ1cZCNDQXZTBfiEspRaKC0qr/FAozKu87fSHQmG+GXMNHwoOjgpXDxPOzmAqzMYh\n/QQOaScwmfNw9XHDq1sHmvXviWfnYNw8nXH1qLtuI1Gz1bEnueXjGH6+exi92zar03PtS9/H6zte\nZ33SenzdfLmj1x1c1fkqTI4XMIlYarwR8LFL4MROY1tAGHSfCqFT0b5dSEs8SlJcbHnYZ51MBkAp\nB1q2DSawYxdade5CYMcu+LYNwcFR/h02JBLs9ai0xGKEf57xIVDxpyjPTKF1e5F1W0F2EYW5Ziy6\n+q/CTs4O1pA3Gd8Gyr4VVNzmaXx7cPW0bnd3wuE8hukJw4b9p7jh/S18e+dQBrZrUS/n3H5yO6/v\neJ0dKTsI8gxiRp8ZTGo/CUeHCwzUjARrn/xiOLbF2ObbBUKnQOhUaNUblCI/O4uTB/dz4kAcyQfi\nOXFwP4U52YDRVx/QvhOBnbrQqpMR9t5+/tKFY0cS7I2AubiUvGMppEdvJHPzH+TEJ2B2cKW0RQCq\nfTd0q2BKXH0otHYXFVo/JLSl+v9mzm5Opz8EPJ1x87R+AHiayh+7eZpw9XC2fiDIh0GZbUfSufrt\nTXx2y2CGd66/uVy01mxI2sDrf7zOvvR9dGrWiRl9ZjA2eCwOygb/bbJPwL6lRsgf+R10KTQLNgI+\ndCoEDQTrcEmtNVknkzlxMJ7k/XGcOBhPyuGDlJrNALh5eePfviMB7Tvi374TAe07Sn99PZJgb4RK\ns7LIWbuWnJWryFu/Hl1cjGPLlniNG4dXRAQegweBkxPFhaUU5pqNsM89Hfjl3xQqPFeQa/xpLqq+\nu8jF2uo3gt8IfDePSh8IHmWPnY3rBw5N73/kXYmZTH3zd96/eQBjQwPq/fwWbSHqSBQLdy7kSPYR\nOjXrxB297yAyJNI2AQ+QlwZxy4yQP7gWLGbwDITQyUbIhwwDxzOHSpaWmDl1NIET++M4efggKYcP\ncupYApZSY7CCi7sH/u064N++Y3noN2/dBocL/dYh/kaCvZGz5OWRu26dcUPUb+vQ+fk4+PjgNXo0\nXpEReAwbhoNL7cdal5hLKcwtoTCvuDzsKwZ/Ya51e97p7aVmS5XHUgpcPE5/A3D1qPSh8LdvCc44\nuzo2+FZdXHIO4xesY9H1/ZgU1spudZRaSllxZAXv7HqHw1mH6ejTsTzgL7iLpqLCLIiPgtifYf8q\nKCkA95bQdZIR8h1GQTXrCZeYzaQdSzCC/sghUg4fIDXhCCXFRYDRjeMX0t7asu+If0gHWgS1xSSL\n0F8QCfYmxFJYSN7GjcbiIWvXYsnOxsHdHc/wUXhFROA5ciQOHrYdnqe1pqTYUiHoiyt9EFgf5525\n3VJa9b8nBwf1ty6h8q4ij8rbjW8GTs4O9fphcORUHuEvR/Pfa3tzRd+gejtvdUotpaxMWMk7u97h\nQOYB2vu05/ZetzOh3QScHGx8A1JxHhxYZfTLx62A4hxw8YYuE4x++U7jwLnmGSgtpaWkH08k5fDB\n8pZ9ypGDFBcUAMYF2matWuMX3A6/4Hb4BrfDL6Qd3r7+cudsLdk02JVS7wPdgV+01rOreN4JOGT9\nAbhXa727uuNJsJ8/bTaTt2WrcUPU6tWUpqWhXFzwGD4c78gIPMPDcfSpZvKpuq5Na8yFpeUhX5Bb\nXP7BUJhrpiCv8jeF4hqvGTiaHM68LlDhG0F13xAcL+B6wfHMAi6Zu4a5V4Zx3aCGs8qWRVtYlbCK\nt3e9zf6M/YR4h3Bzj5uZ2nEqLo510AIuKYJD0cYwyrhfjEnLTO5GuHe/DDpHgmvthoNqi4XMkydI\nPXqEU0ePkJpg/JmZkmwM18SY9Kxl25DTYW/9U26q+jubBbtS6kpgqtZ6mlLqA+AlrfX+Svv0A67V\nWj9am+Ik2G2jfEHvqJXkrFxJSXLy6QW9IyPwGtvwF/TWFk1RQUmFbwDWbqGcyh8G1m3WoafVcXZz\nqhT+Z15Edqv0u4vb6esFablF9J+9iuem9uDmS9rV099A7Vm0hTVH1/De7vf4K+0vWrq25IbuN3BN\n12vO/Uan2iotgYQNRsjvWwq5J8HRGTqMNoZRdp0E7uc+gqi4sIC0Y0dJPXqYU0cTjD8TjlCYl1u+\nj1dLP1oGtaVlUFtatGlLyzbBtAhqi5unly3fYaNiy2B/HVihtV6mlLoOcNNaf1hpnxnA3UAesBu4\nQ2td7f99Euy2py0WY0HvlSv/vqB3ZCReEeOazILepaUWivJKKMgprvID4YxrBTnFNV8vcFDGKCJP\nZ5zdnVh/NJ1uIT707dzSuK/gbyOK6r+LqDKtNduSt/HBng/4/fjvuDu5c3WXq7mh+w3nNlXBubJY\nIHGr9YaoJZB1FJQjtBtuhHy3yeB1/ufXWpObkcaphCPlLfy0pGOkJyWW990DuPs0o2VQ8BmB3zKo\nLe4+zRr8dZwLZctgfx94XWv9p1IqEuintZ5baZ+BQKLW+oRS6hPgO6314kr73A7cDhAcHNw/ISHh\n3N6RqDWtNUVxccbUBpUX9I6MNBb0btvWzlXWL3Nx6eluoJwKF5DzzOUfEPm5xew6mIGvyQkHs24w\nXUQ1iUuP48O/PmTF4RUoFJM6TOLmHjfTpXmXOjlfOa2Nm6DK7npNOwAoaDvYekPUFGNIpS1OZbGQ\nfSqVtKSjpCUeIy3xKOlJx0hLPEZxQX75fq4enrQICqZF6yBatG5D81bGj09AIE6mC7jxqwGxZbC/\nBnyptd5s7ZbpprWeU2kfF611kfXxfYBJa/1KdceUFnv9qnJB79BQvCLG4R0ZiUunTnausGHQWtPx\niWXMCO/EgxFdyruIysL/jAvHeWYKK31A1FUXUW0k5Sbx6d5P+WH/DxSUFDAocBD/DP0n4UHhth1J\nUxWtjekMykL+5B5je6s+1pC/DHxt/2+srIWfnphoDX0j+NOPJ1KQnVW+n1IOePv5lQd981atyx97\n+fo2qmGZtgz2mwB/rfXLSqnngDit9ReV9vkGeBHYA6wE5mitV1V3TAl2+ylOTDJCfmWFBb07dLCu\n9SoLeof+ewU3DAnmyUu7n/NrS0st5cFf8QJxQcXfK3xDqG0XkVsV3UHVjSLKKsrihwM/8NW+rziR\nd4I2nm24rut1XNH5iupXdLK1tIOn73pN2m5s8wstn9qAgB7GmNk6VJiXS+aJ42ScSCLd+qfxcxxz\nYUH5fo4mE80CWhlB37oNzQNb0ywgEJ+AVni1aNngRuvYMti9gfXAamAicB1wtdb6qQr79AS+ABSw\nWGv9ZE3HlGBvGMwpKeSsWkXOypWyoLdVn+ejmNq7Nc9f1rNezmcuLi0P+orDRsuGl5ZtO9dRRK4e\nTmQ5ZHCgMI4E80EszsX0COrG6M6j6NSqXZ13EZXLSjTmk49dDAkbAQ0tOlinNrgM2vSr85CvSGtN\nXmbGGUGfYQ3+zOQT5TddATg6OeHtF2AN+kB8/ANpFtAKn4BAmvkHYnJ1rbe6y9h6uGNzIAJYp7VO\nvtDiJNgbHmNB7zVkR0WRt3ETmM04+fsbd71GRlw0C3oPnrOK8C7+zPu/XvYupUp/G0VU+RtBpfsN\nCnLNNU5nXdddRGfITYF9vxghf3gdWErAO8ga8lMgeAjYsVvEUlpK9qlUsk4mk5WSTGZKMlnJJ4w/\nTyaXz4xZxt2nmRHyAa2soR+Ij18A3n7+eLZoWScTqMkNSuK8lS/oHRVlLOhdWIhj8+bGgt4REU16\nQe+R89fSL7gZC67ra+9SbKasiyg1PZ2VcWvZfGgbuTkFNNMt6ebeg2BTB0zFLmd8GNS2i6hyd1BV\nXUQmlyoCLj8d4lcYXRBXb/IAAA5MSURBVDYHVkNpEXj4GSNrQqdA+5FwITNd1oGC3JzToZ98gqwU\n6+OTJ8k5lYrWp//OlIMDXi398PbzKw97b19/vP0CaN6qNV4tz28uIgl2YRPlC3pHrSQ3OvrMBb0j\nI/AcPrxJLegd8epvdA7wZNH1/e1dSp2xaAsxyTF8F/8dq46uwmwx08evD1d1uYrx7cbj6uhKSbHl\nb91BNXYR5ZqpLkqcTA5/m2LijG8ELqW4ZOzE9Xg0LomrcC1NweTmiuo2yeiX7zAaTPXf7XEuSkvM\nZKemkJ2aSlbqSePxqRSyU0+SlZpCbnpa+Q1ZXQYPY8qsx8/rPBLswuaqXNDbze30gt7hoxr9gt6T\n31hPgJcr708baO9S6kV6YTpLDi7hu/jvOJJ9BE+TJxEhEUzpOIX+Af1rPfnYuXQRlW2vqYvIQZXi\nonJxVdm4OuXj4uWOq18ALgHBuHq74+rhVGFKa2PKalcPE//f3rkGR3Wed/z3gO7yoguS1iJIpsak\n5j5cfCEJFwUkkhS7GTeTxO6XtNPcpnGmmak/pHE+xHEukw/puJ0mLZ0k00kzzsTtOEGubVa4kQ0O\nYCMRLiIWCIwQ4SIEQlpL1mW1Tz6c3ZWAXbG3s7sSz29mZ8++evecv16953/Ovtr3/efn6JpEE4Fx\n/FevMnjlMvlFRdTe9+dJ7ceM3XCV6QO9m7jrow0zMtD7sR+9SUlBHv/9dw9lW0pGUVUOXT7Eb7p+\nQ0t3C8OBYWpLa9lx7w52LN7BvWX3pv2Y4RyDyGqkU1YoHR0KMPLeKKNXLjJytY8R/wijgSJG1ENA\nY9+9R4aKSvOdtLPSKElnpZPbhaHEs5mwSB2YsRsZxAn0/n0k6/XmQG/Ptm3kVVdnW2ZcPL7zAIFg\nkBe+9KFsS8kaw+PD/LbntzSfbmb/xf0ENciK+St4ZPEjbF+0nfnFWVimIhzofWIXEx2vMDI4yCjl\njHg/wkjtJkYr1zASKIqE2UzNMIgn6UzmSOSuP/JcOpmJHM5BDucdFE7JRc7kTGQzdiMrqCojHSci\nE6Iigd5r1zoTonI80PtzP3uLa0Nj7PrKR7ItJSe4MnyFl999mebTzXT2dzJH5vCA9wGaFjWxtX5r\nlkw+dqB3ZGmDslv72EQgOBl1OTTOyHBgyoVg8gJwQ+rZcGDaISNw4i/DJl8YzkMObxeHtksnyz2V\nhZRVT79SZizM2I2so6qMdXUxGBquiQR6r1gRmRCVa4HeX/z5Ic72DbP7a5uyLSXn6LzWye6zu2np\nbuHs4NncMHlVuHRsckLUTYHeLH0UKlPrY8GgMvZ+4JYs5OjbU8sCjL5/67yD+9bVsP3zyc2TMGM3\nco6x7u7IImUjR48CULhkiWPy25so/OAHsz7O+dXnD3P0/HVan2rIqo5cRlU52X8SX7cP31lfxOTX\ne9fTUNfA5rrN1HmytBZRJNB7F1w84pTdvXIyBrDm/ozKCS9lPTLs/LN4dChAQUke1XXJrVBpxm7k\nNOMXL+Jv2YPf52O4rQ1Uyb+n3hmTb2ykaOXKrJj8Uy8cYV9XH/u/vjXjx56JTDX5Pd17ODPgRDIs\nLlvM5rrNbKnbwqqqVe6vVxON/rOhO/nmaQO9ZxJm7MaMIdDXh/+1/8fv8zF08CAEAuTV1kbG5IvX\nrkVcmMUXjad/fYyXj12i/ZuNGTnebKNnsIfW86283vM6bZfbCGiAisIKNi7cyKaFm3jo7ocoLyrP\nvLAbAr33gQZjBnrnMmbsxowkEujta2Fo3z4n0LuqCs/WrU7W64MPIi4uwfrtl07wy7fO0fHMx1w7\nxp3C4Nggv/vj72g938re83sZHBtEEJbOX8qG2g08vOBh1tSscScFajqGrjrJUCd2OUlRwXHw1E7O\neo0S6J0rmLEbM56J94YY2hsr0LuJ0g9/KKFA73j4wavvsPONM3R99xNp3e+dTiAY4HjfcfZf3M+B\nCwc4euUoAQ1QOLeQtTVr2bBgAw/c/QD3V96f/jzX6RgZgJO7nW/YdL12Y6D3sr90ljaIEeidDczY\njVnFtIHeTU3ctXFjWgK9n9tzin/ec5Ku73ycPLdXPryDGRof4tClQxy4eID9F/ZzeuA0AMV5xayq\nXsW6mnWs8a5hVdUqSvKT+2pgwoQDvU/scsx+aqD3skdh8dbbBnq7jRm7MWvRsTEn0LulBf+ePUxc\nu3ZjoHdDA3PnJZcB+u+vn+b7r7zDiWe2U1KQmx/HZyO9w720X26nvbed9svtnOw/iaLMlbksrVzK\nGu8aVlatZPn85dR56tz/x/r4CLz7eloCvdOJGbtxR6ATEwy3tTnfsAkHeufnO4HejducWa+V8Yct\n/3Tfuzzz0gkOf7ORitLZuYLlTMA/5ufIlSMRsz/ed5zRCSf31FPgYdn8ZSyfv9x5VC1nQekC98x+\nYhy630x7oHcymLEbdxzhQO9Bnw+/r4Xxnh4n0Hv9+siEqHyvd9p9/OJgN9948TgH/2kr3nm5vaLg\nncR4cJwz18/QcbWDjr4OOq520NnfSSDozAqdVzCP+8rvcx4VzvOS8iXp/wZOcAJ63pqcEDXQc1Og\n9yPgmb6PpYIZu3FHMxno7WPQ52OsKzSGu3o1nqYmPE2NUQO9/6ftPP/4whHeeKqB+vnZHU81pmds\nYoxT10/R0ddB57VOuq53cer6Kfxj/kidquIqFpcvZtG8RdR56qj31FPnqWOhZyFFeSleuFXhwuFJ\nkw8Hetc/PBkekqZA7zDpTlD6CbAM+D9VfTbZOmHM2I1MM3rmjLNIWUvLDYHe85oa8TQ2RgK9m49c\n4MnnD9PytU0s8SY3O9DIHqpK73AvXde7HKPvP0XX9S7ODZ7DP+6/oW5NSQ31nnoWehbiLfHiLfXi\nLfFSU1KDt8RLeWF5/MM7qtD7h9Cs1+bJQO8FayZjANMQ6J3OzNPHgEdV9XMi8lPge6p6KtE6UzFj\nN7JJJNDb53MCvcEJ9G5qpGPxOv5m3wAvfXUjKz6QofBnw3VUlYHRAXr8PZzzn6PH3xN5nPefp+/9\nPpQbvbBgTgHVJdVUFVdRXlhOWWEZ5YXlt2yX5pdSlFdEcV5x5FF4/TzyTrMzLn+h3dlhzTLH5Jd9\nEryJh6VDeo39X4BXVfVlEfksUKyqP0u0zlTM2I1cYfxyL/7X9uD3tTD8thPo3VdUxlhRcdbXrTEy\niQJBlAmQ0DNBVJxn52dBkGCobjwIqCDOFnMIIihXFxXzqV8cTkplvMYez/e5SoE/hravAWuTqSMi\nXwC+AFBfn95xJ8NIlnxvDZVPPEHlE08Q6O/nyqstdL3oQ8fHsy3NyFGUIEqAIAFUA6HXE6Hn8EVg\n6msNvU8RJhj13uO6xniM/T0gHGp5FxBt1sZt66jqTmAnOHfsCSs1DJfJq6ig9vFP88jjn862FMNI\niXim1rUB4dSB1cDZJOsYhmEYGSCeO/ZfA3tFZAHwceCzIvKsqj49TZ2H0y/VMAzDiIfb3rGr6iCw\nBTgANKjqkZtMPVqdgfRLNQzDMOIhrsUwVLUf+FWqdQzDMAz3seXrDMMwZhlm7IZhGLMMM3bDMIxZ\nhhm7YRjGLCMrqzuKyBWgO4VdVAF9aZKTTkxXYpiuxMlVbaYrMZLVdY+qVt+uUlaMPVVE5FA86yVk\nGtOVGKYrcXJVm+lKDLd12VCMYRjGLMOM3TAMY5YxU419Z7YFxMB0JYbpSpxc1Wa6EsNVXTNyjN0w\nDMOIzUy9YzcMwzBiYMZuGIYxy8hJYxcRr4jsvU2dfBFpFpE3ReRvY5WlWddPRGS/iDw9TZ0vi0hr\n6PF7EfkPEckTkXNTyldmSVtUHSLyLRF5W0T+LUu6ykTkFRHxiciLIlLgZpvFqemWOvG8z01dmW6n\nBHRlvF8loC1b5+O0Hua2f+WcsYtIBfBfOHF70/Ek0KaqHwY+JSKeGGXp0vUYMFdVNwD3isiSaPVU\n9cequkVVtwB7gf8EVgHPh8tV9Vi6dCWiLZoOEVmHE5LyINArItuyoOuvgR+qahNwCfhYNK2Z0hSt\nTgK/i2u6yGA7Jagro/0qEW1ZOh/j8TBX/SvnjB2YAD4DDN6m3hYmlwl+A1gfoyxdTN23j8nEqKiI\nyAcAr6oewgke2SEib4XuMOJaLtkFbdF0bAb+V53/ou8GNmZal6r+SFVbQi+rgd4YWjOlKVqdeN7n\nqq4Mt1PcumJocLNfJaINyPj5GI+HbcFF/8q6sYc+GoU/ErUC/xBnUMfNAdreGGXp0vVkgvv+e+DH\noe23gW2q+iCQD3wiWV0paoumI2faTEQ2ABWqeiCG1nQQz+/rat9KQReQsXZKRJer/SpFbWFcOx9v\nRlUH4/AwV/tYuq9UCaOqX0zyreEA7QGcAO33YpSlRZeIPMftQ73DdecADcA3QkVHVXU0tH0ISOlj\nfAraoumIJ6zcbV2ISCXwr8BfTaM1HSQbzp62dkpBVybbKRFdrvarFLW5fj4miav+lfU79hSIFqDt\nZqh2IvveCBzUyUkCPxeR1SIyF/gkcCSNuhLRFk1H1ttMRAqAF4Cvq2p4cTi32izZcHa3A9tvu/8M\nt1PcumJoyETAfbzHyPT5GA/u9jFVzckH0Dpl+6PAV276+T1AB/AczkerudHK0qhnHk4H+CHwB6AM\nWAY8G6Xud4HHprxeARwFjgHfcaGt4tIWTQfOxf3NUJt1An+WBV1fBvqB1tDjM261WRRNq6Poiab7\nljKX/4bRdGWsnRLUldF+lYi2UL2Mno9TjtMaes64f83omacisgDnCrdbQ2Na0crSeLwKoBF4Q1Uv\npXPfqZKKNhEpBv4CaFfVM7miyy3i0RStjtu/Sy62FSSvy81+laq2XMBN/5rRxm4YhmHcykweYzcM\nwzCiYMZuGIYxyzBjNwzDmGWYsRuGYcwyzNgNwzBmGX8Cl4HbYn0q8kkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1bff865d940>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签\n",
    "plt.rcParams['axes.unicode_minus']=False #用来正常显示负号\n",
    "\n",
    "#构造损失函数\n",
    "def e1(s):\n",
    "    if s>0:\n",
    "        return 0\n",
    "    else:\n",
    "        return 1\n",
    "\n",
    "def e2(s):\n",
    "    return max(0,1-s)\n",
    "\n",
    "def e3(s):\n",
    "    t=max(0,1-s)\n",
    "    return t**2\n",
    "\n",
    "def e4(s):\n",
    "    return max(0,-s)\n",
    "\n",
    "def e5(s):\n",
    "    return 1/(1+np.exp(s))\n",
    "\n",
    "def e6(s):\n",
    "    return np.exp(-s)\n",
    "\n",
    "x=np.arange(-1,1,0.01)\n",
    "\n",
    "y1=[e1(i) for i in x]\n",
    "y2=[e2(i) for i in x]\n",
    "y3=[e3(i) for i in x]\n",
    "y4=[e4(i) for i in x]\n",
    "y5=[e5(i) for i in x]\n",
    "y6=[e6(i) for i in x]\n",
    "\n",
    "plt.plot(x,y1,label='e1')\n",
    "plt.plot(x,y2,label='e2')\n",
    "plt.plot(x,y3,label='e3')\n",
    "plt.plot(x,y4,label='e4')\n",
    "plt.plot(x,y5,label='e5')\n",
    "plt.plot(x,y6,label='e6')\n",
    "plt.legend()\n",
    "plt.title('损失函数比较')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "因此e2,e3,e6,即(a),(b),(e)为$[\\![y\\neq sign(w^Tx)]\\!]$的上界"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 4\n",
    "\n",
    "这题要找哪些函数不是处处可导，(d),(e)显然是处处可导的，因为$f(s)=max\\{0,s\\}$在原点不可导，所以(a),(c)不可导，唯一有点疑问的是(b)，我们来看下$f(s)=(max(0,s))^2$\n",
    "$$\n",
    "f(s)=\\begin{cases}\n",
    "s^2(s\\ge0)\n",
    "\\\\0(s<0)\n",
    "\\end{cases}\n",
    "$$\n",
    "显然这个函数也是可导的，所以由复合函数的性质可知(b)也可导。因此这题选(a)(c)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 5\n",
    "\n",
    "$err(w) = max(0, -y_nw^Tx_n)$的意思是对于分类正确的点$err(w)=0$，对于分类不正确的点$err(w)=-y_nw^Tx_n $，我们来求梯度(不考虑不可导点)\n",
    "$$\n",
    "\\frac {\\partial(-y_nw^Tx_n)}{\\partial w_i}=-y_nx^{(i)}_n(x^{(i)}_n表示x_n的第i个分量)\n",
    "\\\\\\nabla(-y_nw^Tx_n)=-y_nx_n\n",
    "$$\n",
    "所以对于分类错误的点$(x_n,y_n)$，根据SGD,更新规则为\n",
    "$$\n",
    "w(t+1)=w(t)+\\eta(-\\nabla(-y_nw^Tx_n))=w(t)+\\eta y_nx_n\n",
    "$$\n",
    "所以PLA可以被看成$err(w) = max(0, -y_nw^Tx_n)$的SGD"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 6\n",
    "\n",
    "$$\n",
    "\\frac{\\partial E(u,v)}{\\partial u}|_{(u,v)=(0,0)}=e^u+ve^{uv}+2u-2v-3 |_{(u,v)=(0,0)}=-2\n",
    "\\\\ \\frac{\\partial E(u,v)}{\\partial v}|_{(u,v)=(0,0)}=2e^{2v}+ue^{uv}-2u+4v-2|_{(u,v)=(0,0)}=0\n",
    "$$\n",
    "\n",
    "所以\n",
    "$$\n",
    "\\nabla E(u, v)=(-2,0)^T\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 7\n",
    "\n",
    "编程处理即可,需要利用刚刚的偏导数公式，答案为2.82500035668"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2.82500035668\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "def E(u,v):\n",
    "    return np.exp(u)+np.exp(2*v)+np.exp(u*v)+u*u-2*u*v+2*(v*v)-3*u-2*v\n",
    "\n",
    "def partial(point):\n",
    "    u=point[0]\n",
    "    v=point[1]\n",
    "    pu=np.exp(u)+v*np.exp(u*v)+2*u-2*v-3\n",
    "    pv=2*np.exp(2*v)+u*np.exp(u*v)-2*u+4*v-2\n",
    "    return np.array([pu,pv])\n",
    "\n",
    "point=np.zeros(2)\n",
    "eta=0.01\n",
    "\n",
    "for i in range(5):\n",
    "    point-=eta*partial(point)\n",
    "    \n",
    "print(E(point[0],point[1]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 8\n",
    "\n",
    "这题需要用到多元泰勒公式，可以参考[维基百科](https://zh.wikipedia.org/wiki/%E6%B3%B0%E5%8B%92%E5%85%AC%E5%BC%8F#%E5%A4%9A%E5%85%83%E6%B3%B0%E5%8B%92%E5%85%AC%E5%BC%8F)\n",
    "$$\n",
    "b_{uu} =\\frac 1 2\\frac{\\partial^2 E(u,v)}{\\partial u^2}=\\frac 1 2\\frac{\\partial}{\\partial u}\\frac{\\partial E(u,v)}{\\partial u}=\\frac 1 2\\frac{\\partial (e^u+ve^{uv}+2u-2v-3) } {\\partial u}=\\frac 1 2(e^u+v^2e^{uv}+2)\n",
    "\\\\b_{vv}=\\frac 1 2\\frac{\\partial^2 E(u,v)}{\\partial v^2}=\\frac 1 2\\frac{\\partial}{\\partial v}\\frac{\\partial E(u,v)}{\\partial v}=\\frac 1 2\\frac{\\partial{(2e^{2v}+ue^{uv}-2u+4v-2)}}{\\partial v}=\\frac 12 (4e^{2v}+u^2e^{uv}+4)\n",
    "\\\\ b_{uv}=\\frac{\\partial^2 E(u,v)}{\\partial u \\partial v}=\\frac{\\partial(e^u+ve^{uv}+2u-2v-3)}{\\partial v}=(e^{uv}+uve^{uv}-2)\n",
    "\\\\b_{u} =\\frac{\\partial E(u,v)}{\\partial u}=e^u+ve^{uv}+2u-2v-3\n",
    "\\\\b_{v}=\\frac{\\partial E(u,v)}{\\partial v}=2e^{2v}+ue^{uv}-2u+4v-2\n",
    "\\\\b=E(u,v)\n",
    "$$\n",
    "将$u=v=0$带入可得\n",
    "$$\n",
    "b_{uu}=\\frac 3 2\n",
    "\\\\b_{vv}=4\n",
    "\\\\b_{uv}=-1\n",
    "\\\\b_u=-2\n",
    "\\\\b_v=0\n",
    "\\\\b=3\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 9\n",
    "\n",
    "由题设我们知道Hessian矩阵正定，此处的Hessian矩阵为\n",
    "$$\n",
    "\\nabla^2 E(u, v)=\n",
    "\\left[\n",
    " \\begin{matrix}\n",
    " \\frac{\\partial^2 E(u,v)}{\\partial u^2} &   \\frac{\\partial^2 E(u,v)}{\\partial u \\partial v}  \\\\\n",
    "   \\frac{\\partial^2 E(u,v)}{\\partial u \\partial v} &  \\frac{\\partial^2 E(u,v)}{\\partial v^2}\n",
    "  \\end{matrix}\n",
    "  \\right]\n",
    "$$\n",
    "那么最优方向为\n",
    "$$\n",
    "\\left[\\begin{matrix} \n",
    "\\Delta u^*\n",
    "\\\\ \\Delta v^*\n",
    "\\end{matrix}\\right] \n",
    "=- (\\nabla^2 E(u, v))^{-1}\\nabla E(u, v)\n",
    "$$\n",
    "一般性的结论可以参考凸优化的课本，关于这题可以简单证明下。\n",
    "\n",
    "现在要对$\\hat E_2(\\Delta u, \\Delta v)=b_{uu}(\\Delta u)^2+b_{vv}(\\Delta v)^2+b_{uv}(\\Delta u)(\\Delta v)+b_u\\Delta u + b_v\\Delta v + b$求最小值，令$\\Delta u =t,\\Delta v=s$\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\hat E_2(\\Delta u, \\Delta v)&=b_{uu}(\\Delta u)^2+b_{vv}(\\Delta v)^2+b_{uv}(\\Delta u)(\\Delta v)+b_u\\Delta u + b_v\\Delta v + b\n",
    "\\\\&=b_{uu} t^2+b_{vv}(s^2)+b_{uv}st+b_ut+b_vs+b\n",
    "\\end{aligned}\n",
    "$$\n",
    "由Hessian矩阵的正定性我们知道\n",
    "$$\n",
    "2b_{uu}= \\frac{\\partial^2 E(u,v)}{\\partial u^2}>0\n",
    "\\\\ 2b_{vv}= \\frac{\\partial^2 E(u,v)}{\\partial v^2}>0\n",
    "\\\\ \\frac{\\partial^2 E(u,v)}{\\partial u^2} \\frac{\\partial^2 E(u,v)}{\\partial u^2}-  (\\frac{\\partial^2 E(u,v)}{\\partial u \\partial v})^2=4b_{uu}b_{uv}-b_{uv}^2>0\n",
    "$$\n",
    "所以可以对$\\hat E_2(\\Delta u, \\Delta v)$进行配方(二次项系数不为$0$)，配方得\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\hat E_2(\\Delta u, \\Delta v)&=b_{uu} t^2+b_{vv}(s^2)+b_{uv}st+b_ut+b_vs+b\n",
    "\\\\&=b_{uu}(t-a)^2+b_{vv}(s-b)^2+b_{uv}(s-a)(t-b)+C\n",
    "\\end{aligned}\n",
    "$$\n",
    "其中$a,b,C$均为常数，后续会求解出来，令$t_1=t-a,s_1=s-b$\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\hat E_2(\\Delta u, \\Delta v)&=b_{uu}t_1^2+b_{vv}s_1^2+b_{uv}t_1s_1+C\n",
    "\\\\&=b_{uu}(t_1+\\frac{b_{uv}}{2b_{uu}}s_1)^2+(b_{vv}-\\frac{b_{uv}^2}{4b_{uu}})s_1^2+C\n",
    "\\end{aligned}\n",
    "$$\n",
    "之前已经有$b_{uu}>0,b_{vv}>0,4b_{uu}b_{uv}-b_{uv}^2>0$，所以\n",
    "$$\n",
    "b_{vv}-\\frac{b_{uv}^2}{4b_{uu}}=\\frac{4b_{uu}b_{vv}-b_{uv}^2}{4b_{uu}}>0\n",
    "$$\n",
    "从而\n",
    "$$\n",
    "\\hat E_2(\\Delta u, \\Delta v)\\ge C\n",
    "\\\\当且仅当t_1+\\frac{b_{uv}}{2b_{uu}}s_1=0,s_1=0时等号成立\n",
    "\\\\即t_1=s_1=0\n",
    "$$\n",
    "因为$t_1=t-a,s_1=s-b,\\Delta u =t,\\Delta v=s$，所以等号成立条件为\n",
    "$$\n",
    "\\Delta u= t=a,\\Delta v=s=b\n",
    "$$\n",
    "接下来求解$a,b$\n",
    "$$\n",
    "\\begin{aligned}\n",
    "b_{uu}(t-a)^2+b_{vv}(s-b)^2+b_{uv}(t-a)(s-b)+C&=b_{uu}(t^2-2at+a^2)+b_{vv}(s^2-2sb+b^2)+b_{uv}(st-as-bt+ab)+C\n",
    "\\\\&=b_{uu}t^2+b_{vv}s^2+b_{uv}st-(2ab_{uu}+bb_{uv})t-(2bb_{vv}+ab_{uv})s+C^{'}\n",
    "\\\\&=b_{uu} t^2+b_{vv}s^2+b_{uv}st+b_ut+b_vs+b\n",
    "\\end{aligned}\n",
    "$$\n",
    "那么\n",
    "$$\n",
    "\\begin{cases}\n",
    "-2b_{uu}a-b_{uv}b=b_u\n",
    "\\\\-b_{uv}a-2b_{vv}b=b_v\n",
    "\\end{cases}\n",
    "\\\\\n",
    "\\left[\n",
    "\\begin{matrix}\n",
    "-2b_{uu} &   -b_{uv}  \\\\\n",
    "-b_{uv}  &   -2b_{vv} \n",
    "\\end{matrix}\n",
    "\\right]\n",
    "\\left[\\begin{matrix} \n",
    "a\n",
    "\\\\ b\n",
    "\\end{matrix}\\right]\n",
    "=\n",
    "\\left[\\begin{matrix} \n",
    "b_u\\\\b_v\n",
    "\\end{matrix}\\right]\n",
    "\\\\\n",
    "$$\n",
    "回顾之前的等式\n",
    "$$\n",
    "2b_{uu}= \\frac{\\partial^2 E(u,v)}{\\partial u^2}, 2b_{vv}= \\frac{\\partial^2 E(u,v)}{\\partial v^2},b_{uv}=\\frac{\\partial^2 E(u,v)}{\\partial u \\partial v}\n",
    "\\\\b_{u} =\\frac{\\partial E(u,v)}{\\partial u},b_{v}=\\frac{\\partial E(u,v)}{\\partial v}\n",
    "\\\\\n",
    "\\nabla^2 E(u, v)=\n",
    "\\left[\n",
    " \\begin{matrix}\n",
    " \\frac{\\partial^2 E(u,v)}{\\partial u^2} &   \\frac{\\partial^2 E(u,v)}{\\partial u \\partial v}  \\\\\n",
    "   \\frac{\\partial^2 E(u,v)}{\\partial u \\partial v} &  \\frac{\\partial^2 E(u,v)}{\\partial v^2}\n",
    "  \\end{matrix}\n",
    "  \\right]\n",
    "  \\\\\n",
    "  \\nabla E(u, v)=\n",
    "\\left[\n",
    " \\begin{matrix}\n",
    " \\frac{\\partial E(u,v)}{\\partial u}  \\\\\n",
    "\\frac{\\partial E(u,v)}{\\partial v}\n",
    "  \\end{matrix}\n",
    "  \\right]\n",
    "$$\n",
    "原方程可化为\n",
    "$$\n",
    "-\\left[\n",
    "\\begin{matrix}\n",
    "\\frac{\\partial^2 E(u,v)}{\\partial u^2}&   \\frac{\\partial^2 E(u,v)}{\\partial u \\partial v} \\\\\n",
    "\\frac{\\partial^2 E(u,v)}{\\partial u \\partial v}  &    \\frac{\\partial^2 E(u,v)}{\\partial v^2}\n",
    "\\end{matrix}\n",
    "\\right]\n",
    "\\left[\\begin{matrix} \n",
    "a \\\\\n",
    " b\n",
    "\\end{matrix}\\right]\n",
    "=\n",
    "\\left[\\begin{matrix} \n",
    "\\frac{\\partial E(u,v)}{\\partial u} \\\\\\frac{\\partial E(u,v)}{\\partial v}\n",
    "\\end{matrix}\\right]\n",
    "\\\\   \n",
    "\\left[\\begin{matrix} \n",
    "\\Delta u^*\n",
    "\\\\ \\Delta v^*\n",
    "\\end{matrix}\\right] =\n",
    "\\left[\\begin{matrix} \n",
    "a \\\\ b\n",
    "\\end{matrix}\\right]\n",
    "=- (\\nabla^2 E(u, v))^{-1}\\nabla E(u, v)\n",
    "$$\n",
    "这样就验证了牛顿方法的正确性。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "#### Problem 10\n",
    "\n",
    "有了公式之后编程实现即可\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2.36082334564\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from numpy.linalg import inv\n",
    "\n",
    "def E(u,v):\n",
    "    return np.exp(u)+np.exp(2*v)+np.exp(u*v)+u*u-2*u*v+2*(v*v)-3*u-2*v\n",
    "\n",
    "def partial(point):\n",
    "    u=point[0]\n",
    "    v=point[1]\n",
    "    pu=np.exp(u)+v*np.exp(u*v)+2*u-2*v-3\n",
    "    pv=2*np.exp(2*v)+u*np.exp(u*v)-2*u+4*v-2\n",
    "    return np.array([pu,pv])\n",
    "\n",
    "def dpartial(point):\n",
    "    u=point[0]\n",
    "    v=point[1]\n",
    "    puu=np.exp(u)+np.exp(u*v)*(v**2)+2\n",
    "    pvv=4*np.exp(2*v)+np.exp(u*v)*(u**2)+4\n",
    "    puv=np.exp(u*v)*(1+u*v)-2\n",
    "    return np.array([[puu,puv],[puv,pvv]])\n",
    "\n",
    "point=np.zeros(2)\n",
    "eta=0.01\n",
    "\n",
    "for i in range(5):\n",
    "    point-=inv(dpartial(point)).dot(partial(point))\n",
    "    \n",
    "print(E(point[0],point[1]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 11\n",
    "\n",
    "这题问的是这六个点能否被直线和二次曲线shatter，比较简单，作图即可，这六个点可以被shatter."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "#### Problem 12\n",
    "\n",
    "先回顾题意，假设有$N$个点$x_1,...x_N,x\\in R^d$，现在构造这样一个$R^d$到$R^N$的映射\n",
    "$$\n",
    "(Φ(x))_n = z_n = [\\![x = x_n]\\!]\n",
    "$$\n",
    "这题有点抽象，我们举$N=3$的例子看一下\n",
    "$$\n",
    "\\begin{aligned}\n",
    "Φ(x_1)&=((Φ(x_1))_1,(Φ(x_1))_2,(Φ(x_1))_3)\n",
    "\\\\&=([\\![x_1 = x_1]\\!],[\\![x_1 = x_2]\\!],[\\![x_1 = x_3]\\!])\n",
    "\\\\&=(1,0,0)\n",
    "\\end{aligned}\n",
    "$$\n",
    "所以其实这题很简单，是将$N$个$R^d$空间上的点映射到$R^N$上，我们知道$R^N$上感知机的$d_{vc}=N+1$，所以$N$个点一定能被shatter,所以\n",
    "$$\n",
    "d_{vc}(H_Φ)=+\\infty\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "这里先根据题意产生一组点，作图看一下"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD6CAYAAACoCZCsAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJztnX+UXkWZ5z9lQu/SzY8hMSeObgLC\nJGhGh6z2xqgoZAieIcvMchgSXM0MMOshg47LqLMHURwXDwkyc8hO64xAHDKAHJ2TzLqO4I9pyS8R\nCNDMiSCRJiGiEZCEblTSPaYnsfaP977p27fvj7p1q+rWfd/6nnNP377vvVXPrVv1PFXP89TzCCkl\nAQEBAQHdh1fVTUBAQEBAQD0IAiAgICCgSxEEQEBAQECXIgiAgICAgC5FEAABAQEBXYogAAICAgK6\nFEEABAQEBHQpggAICAgI6FIEARAQEBDQpZhZNwF5ePWrXy1PO+20uskICAgIaBQee+yxl6SUc4ru\n81oAnHbaaQwNDdVNRkBAQECjIIT4scp9QQUUEBAQ0KUIAiAgICCgSxEEQEBAQECXIgiAgICAgC5F\nEAABAQEBXYogAAICAgK6FNoCQAgxVwhxf87vxwkh7hFCPCCE+JOsawEBAQEB9UBLAAghTgHuBPpy\nbvsw8JiU8p3AJUKIEzOuBQQEBATUAN0VwFHgUuCXOfecC2yKzr8L9GdcmwIhxJVCiCEhxNDBgwc1\nyauG0bEJbtvxDKNjE7XUHxDQyfBpfOXRYpXOsRF4YKD111WdKdASAFLKX0opf1FwWx/wXHQ+CszN\nuJYse4OUsl9K2T9nTuFOZivYPLSfG7/1FJuH9tdSfzfAJybgCi7f2Yf2zaKhPb7ufPBHzmlM0pQ3\n1q3ygV13w3f+svXXVZ0psBkK4hBwPPAL4ITo/7Rr3mFl/7wpf0fHJtg8tJ/li+Zy3+4XWdk/j1l9\nPcbqa5dvqlzT5dmor93RAdacc4bTuuuCqXf2ra6yNLTH1fjE0co0lv3+SZqSYz2OvN8qY/HqqX9d\n1JkCmwLgMeBs4J+As4CdGde8w6y+nikdst1pdu4bYdtwSy1lclCZHqyuB79OffGOXoWJ23hXJXrG\nRlqzt8WroW+2UrkuB7epuqp8mywa2uNrdGyC3p4ZlWgs+/2TNCXHehqdVtA3G955tds60yCl1D6A\n7dHf3wX+LPHbqcCTwADwKDAj7Vpe+W9961ulDxg5dFjeun2v3HvgFXnr9r1y5NBhK+WbKtd0ebbr\nu3X7XnnqNffKW7fvdV63Nj3f+xspP31S628Ho8q3qQLV7+q6rzcFwJBU4eEqN+kewGuBVcDJedey\nDl8EQIBdOBnEh15qMetDLxXWqURPrLxOZkJp7+bifduC5/KND+t/owaj6vupCgCrG8GklM9LKTfJ\nmME47Vo3o6yxropxz6Rh0GRZs/p6WNk/j81D++0ZBBNGtzxjW3sZnqvyaC/h+2Z75TRg2vib1hYu\n3ndl/zyWnTmHbcMHnRlo6zSclzFOm4TX+QC6AWV1mKn3K+qjTerLG2e3WLwaJsZhYgzGRozq410b\n7vLgwv6T9b5Z9gIdO8Ksvh5uXrX42HN5NJhyBKjTcF7GOG0UKsuEug4XKqC9B16Rl298WO498Ir1\nukzUn7o0jPTRY1tvzl02mlw2N9JuYVBv76MKYuTQYbl+8Cm5fnC4Frqy7AV5dgQT7ahip1Cpp85v\napoPoagC6voVwA337o48e3bzD1cscV7/fbtfZNvwQZae/iJnnHNC4f2pXgKRK9nm8Xdy45bsGYyq\nh4HKjMq0t4IT74cM17s0FLWB7dmirmvtwJa9XHvBGwqfseE+mzVrzZvNmmhHldmySj2m+qBO25bl\nA6bQ9QLgugsXAbujv+5hZKkX6aN/f2yCX/WkL5nL4M4Hn2Vgyx5GxiaYHennqzIJL/z1M1zv0lDE\nMGwv0afVr6DmK0OTDQGWxUDzGKuJdlRh3C7VdFXdop1CZZlQ12FCBeTjUl0JkZfJ6IHnndO/fvAp\neeo198r3bXjImAtgXe6Euqi730yr37Dbad3v5wp1vKcPbUtQAbXgw45ILUReK8P7RrjxybcB7ui/\n7B2vp7dn5pSdz1Xhk6FUBc435BTVX0J9pVW+BnzbcZ72ex3jX7dt61gld7wAaBrjOYZooJ+5cBXX\nzh93Sn+8A5vSR2YOCo0dtVXhhTqqLEqor1zBtx3nab/7Ov59EVa1q3nyDucbwRKbhQIcoIYdtWXU\nUT4s503CS08wxXFXVF8ZeuL31vGN0/qgSTrwYSOYz0jd9JERoa90OR7B+wiUi1fD+Z+ZotqwTfPK\n/nlce8EbShlLfdjkZQL3PPg4I4N/zT0PPl6qndPuVdowp1Ke4ribxSusmXkvs3gl/fcS9MS/ax3f\nOK0PVmlPXXS8CigLqcstDT2r7zYG7yNQpqg2dGlWVe2U0dHWpULIe5dKAdpmfpfe477C+MyFfGno\nFOV2trrxr19x3LUFBVRWh6V917rUrLVCZZlQ12FTBWRqueV8+VhSTeWSPlttqlpukWqnSQHG8t6l\n0ntqxjDyYuNfTSpaH/pDWeBDMLiqh0sbQGM+cpdEoYxDVWdf9A1Vy/HBZbVKsDpT7RXQwvrBYXnq\nNffK9YPD036rHFTQElQFQNeqgJLwXZVzDLrugDV42yijgDZVNUzRslo1fowPniNV4tSr0t+YPl87\nZOLvJPLasAnt2zUCoEhv6sOgV4KuO2CO/rR2t8gC3a4pfWm8nNt2PJM5OL3Rz2pClf7G9Pma0d4X\nkxb8bnziKFe+6/WMTxxldGxiyvhpQvt2jRdQkaXflAXeW6+gFG+bNsp4QVT1HClLmy20vTCWL5rr\n5/eKQaUddfpdHV4nPiOrDbPaqRV7aQ97DhxiYMueaeMn7Tnf+EPXrABcSeOiZZ/t/MJZ9a3sn8es\njJWDrRgyyvcmVjUuViTtwZm3EvAFKu3YBHUDeLDazIFuesnli+ay9HS1HfPx9LI3r1pcextoCwAh\nxO3AIuAbUsobUn6/Crg0+vc3gIeBDwH7ogPgw1LKJ3RpKIMqy3rVTttaEh7h6vMWZHYG2/mFs+oz\nVU8ZYaErdF0yMy+X6QmbiAqNXr5HCnwWVGXbUGfH/Mr+ecfG/uah/fW3gYqlOHkAFwN3ROcbgQUF\n938e6AfeAtykWo8vKSHLeo+sH3yq0DPARH5hUzHOffB6iaNW7xRTroZVyulgT6+m5au2ARc0YdMN\nFPgcsCI6fy9wRc69rwM2RecfpJUU/hHgdmBmyv1XAkPA0Pz58601UBmU8R9fP/iUfP8XdzphqKYY\nt4+DpDaoMF8V5l6WicfLbFhIEt19Gy7genJTd4iJNlQFgK4KqA94LjofjWb2WfgQcEt0/iiwXEr5\nghDiLmAF8PX4zVLKDcAGgP7+/ul+VxaRpepRVR/N6uuht2cm39v7EsvOnGN9OZ5csibpL7sztm2g\n8lE/6wwqbrYqO1LLuusmy7QY+M20Hv7OB3/EwJa9jE8c4SPnn1lJzfPywRd44t4v8OTcC7n0nP9c\nmT7XqrH4uwPeqrva0BUAh4Djo/MTyPAmEkK8ClgGfDK69LiU8nB0PgQs0KzfCkxnJ7LNRJOCKUm/\nkXzDPsLmngYVN1sV5l7WXdfh/g7z31lM+VuF6Q5/+xbe/ePP8b29B9nc9+rK9Ll26a07xERpqCwT\nkgfwx8BfROfXA+/LuO8c4G9i/28CzgJmAFtprQa8sQH4tHTVQdWluIv3r1zHoZekvPuSjtWRl4aG\nvcCLsA4ZGD3wvNyx8Tr5hW/sbOw49AEoqoBE695yEEKcBNwPbAEuoGUHWCmlvC5x37qIkK9G/78J\n+DKtqcLXpZSfJAf9/f1yaGioNH3dBJ/d6tLQdru89oI36M3MHhhoqUoWvAcuutW/Xc2u4fMO74po\nWt/2CUKIx6SU/UX3aamApJS/FEKcC5wP/JWU8mfA91Pu+0Ti/x8Av6NTZxV0ckfKW877+N6VdbJx\nVUmHMTwt6O4Mb4DgaIxKsipq/BbaO4GllC9LKTdFzN9rdFpM9zjydrT6+N6Vd5+2GZ6nTKsx0Mh9\n4Rpl8jZkwcjO27GR1spzbKRyeabykJhCV+wE9nWTjIkZet6OVl/fu6vg60zbcI5hG+9pwoBrZOdt\nzENr85ELK61KTOUhMQYVQ0Fdhy8bwWxBZeOYKrxM9ecDdHzqk880YVOXjb0DZcqs8J57D7wiL9/4\nsNx74BUNIvMxcuiwvHzjw9X2AmjmUMiix8WeCUJKyGwkl2F1BWhqL3FBVFbVmAzs5aPqSBs6y+vk\nM1WW6AtXtAzWC1dM/y2mWsiEyj1VaTRRZoWAfjfcu5ttwwe54d7dpZ8tGruz+nq4edXiaqqkmNqx\n6jhLPl/3WOsKFVASVf3l26iqwolvwAIYnzgyLaRsHVjZP4/xiSOpIW6PwVfVRhI6y+vkM1WW6E9/\nE/YMwmlnw5yEsVZlQ5lqGkQbaoQyZeoao4HrLlwE7I7+loPK2LW5F6AqD6hdTauyTKjrsKUCMrUM\nM7nN3Ld4PIX0mFZtNCz8gTLy3kvlnU2osDoYNtWVTYyj1QYhJWQxTOvz6qTFNArpMc1kOjgAmhJM\ntqdHbelbvy4DFebu6/upCoCuVAG1YcLPeHziCHc++CyXveO0Sqob37JQFdJTYcmfijo8IXxSYyVV\nPVVoM9mWFduoyb78KuoZ38ZtWXS1ANDVv7X1fuMTRxnYsheA3p4ZVjuCj5u6jMKkQFFlWqr6dRdI\nMu0qtJlsy4ptVLuOOwVlgySWKffOB38EiMoTQlfoSi+gNmb19bCyfx6bh/aX8gCanNVIrj7vt3IT\nwJRBnkdD3d4Cyt4oRc9WKUcVqt4rSc8VF7RlIbnBrcirxhWtFdN1qnrNuPTMszWWWiki96amh/QV\nHb0CiEt6IFXql1mixtM5gn7Ez6wZSB4ttc+kqswE48+C/Vm3qgokOVP2aUVQNIt3RWseHQZVaKY8\n81Rgayy1vedAeLXiyUNHCwCV2Ny28uG2kcbss8rJo2UWr7Bm5r3AaqAGfXUZvXKSMaQ9a1PXr6sC\nqXNHZllUodUU4zYohJJ9vyyTLqMitaW3n9XXw0fOP9N4uVahYimu66jqBWQ6O49OGWmeBFq0RJ4d\nY1tv9tLrYAo88kIJSIGp7+ORu6lv7ph1ewcRvICmS3rXySVaSeKPcvV5vzVlJqM1A4lmepvH38mN\nW1qrh7b9InfWU4enS5XVQifDl3c1tdIx7QlWAbWrSBNoivdTRwuAutEyCu3h2gveUN0jIBpsvz82\nwa969h9j/oWdrA69tgpjaDPDiXHY8dnWtar0+cJgs+CLjcEjxm0Kvrlj+iaQstDVXkAqngZVvBFM\nhLNNIu5VoVT+4tVwzsdhYqwe75YsHGOGspKXSWqZdYU4LvLMqehRUwm6XkNZz1n2QrIdn8t2+boe\nhq6hLQCEELcLIR4SQlyX8ftMIcRPhBDbo+PN0fXrhRCPCiH+TrduU1BxB6viMmYyQJt2+X2zoacX\ndtykzxhtDPY2M1yyplx8/zxa6mSwUCyA6sxloCscs57LuG6KsSbHnWmGnTeubb2Dl1AxFCQP4GLg\njuh8I7Ag5Z63ADclrr2VVhpJAXyamnMCqxhq6jbmGEFVY52q0TBejy0Doc8G5rqMojoxhVRpzbov\n47opY2xy3Jk28uaN6/WDw1GY9mFjdbjmI9iMBQR8DlgRnb8XuCLlng8CTwKPALfTsjd8FPhg9PtS\n4Pq8ejo9H4BXKApatnWtlFvX5TOMOHNOMmpTzNEjzxNvoNPWGYLUp/hYtsvNKnP94FPH8nSYqt+1\nl5KqANA1AvcBz0Xno9FsP4lHoxn+C0KIu4AV0XPPxJ6bm3xICHElcCXA/PnzNckzg44PvxBH3EC5\neHXr/4UrWuGMF6+Gnr7W7z292QbEPH9/UwbQDjRgVoZOGIkMT6Cq3is2fexNl5v1rpe94/X09syc\nYlur2i6+GoV1BcAh4Pjo/ATSbQmPSykPR+dDwAKV56SUG4ANAP39/VKTPiNoiisXUN0DJs4Q2gzk\n2e+1Ytknf89CkjnHz5u0ySoLNr2MqpSdbHedbxXBV0ZlA1nvGhc2abv/deCbl9IxqCwTkgfwx8Bf\nROfXA+9LuWcTcBYwA9gKLAfeDfxt9PsVwCfy6qlbBWRj2ZlMf6daR+F9JnXjbRXCgadL6X+9Q5NC\nVvts2+hiTFHdNKXfS/sqoK8B9wshXgtcALxXCHGDlDLuEfQZ4Mu0DL5fl1LeJ4R4FXCjEGIA+L3o\n8BY2pHY7/R3s5h+uWKK8yii8z2RogPjsMJnFCvzxZy+CaTptrmI6YYXUgZiySth1i3p/8n1PShsq\nUiLtAE4BVgGvKfnc8cAlwOlF99axArBtrbe2AqiCretas8+t69J/N5kk3SV8pdOFt1SAeZT5VjWv\n6LCdFF5K+bKUcpOU8mcln/s3KeU/SSn36dZtE+2Z9sc27bKygeOMOSfwD1cs4Yw5JwDqewXs7imQ\nib8JJH2+q/izp/nx29pUVKfffR7i7Vn35rWmQKGP5PnvG/HtT/Sn3DLr3pOiiBAKIkLc2LNz3wjb\nhg+yeWi/kgqocd5CyeXpkjUtL5+szmpSPZGmlmmKSskUXEdH7QQo9JE8NakNh47cMpviraayTKjr\ncLkRLG7sKatu8S0SoZQyf7mqsjx1mbDcohokbOSruU6H5eR9a5f7CHwAISl8Maow/Thsd4Qq4aNT\nmXyFjUKl77GBEkzFS+FcFnW0s6k6XdEebClToCoAukYFlKamiVv4q3j82Pbx1Vq+tlUKC1e0dKdx\nbwSV5amK2qcuz5USKiMjfu11eHSMjcAjGwAJb1rZula1ncu8h6lv66qPJPtEU7xw6oaKlKjrMLkC\naPJMsNIKoxP9y3Vne7rPlW1DE7PRdp0mv13ZmE5Z+0B0YXOWnizbQb/vBBVQ16wAfNzhqGo8rrTC\nsD0Dq2OmpWpgS9Kma2wu24YmjNqLV7dyJSDNfbuFK1q7uxeuyL8vbSe4CYOmTWO/zm7oimhUpIAs\nqEiJuo66dgK7kuxNXpUcg85My5W+1lZAugiZ/cRXfXTRt0rO/Ju0AjCMpkcKJhiB9WErpG3Z35VR\n58DSqbvq8rxqKGNDaJwAL2qPMt+lQcxcB6rfVncM2xYeqgKga1RAcRSpXkypi4qWiMaMx3X60ev4\nO+ssz+PqHNX3LUtbSXWWj2rFXBS1R1vt1M4el9cGKt8grz1V27omY67qt9VRA42OTfCxTbuikDD1\nqo+6UgC4YszOGESex08SPnhH6AiNZLjq+F9TKClIvY3wqIt29rjv/GVrY2BJYTFtYpXXnqptHbdH\nXHSrkT6rYntT/bY6Y3zz0H62DR9k2Zlz6p88qCwT6jpsqYB81t1VQlXffd1NWi7UAZ1SRxX41gaJ\nvjRNbVKUZCivnrg94u5LjHr01K26c8F/CDaALkRVZp0mHHzeENYpUGW6ptu5qkBJPF/I2HSFiyvj\nfQchCACHcNqhTA4GlSifCvWNHnhePnTXp+Togeer09SNKOufXzX0Rtl6TSEYmZ1BVQBoRwMNmETb\nprB5aL/5wpNREE1Gj1SJ8qkQUXPT7nHe++Tb2LR7vDpNbdiKEFon2u90cM/Ud1ONHKkT3TSvv7iO\nWNmur22ryvu2vkZyLQEjEUgtoyuNwHnQiexp1dibNJaZNIDmlVXCWGzl/TsxQmjWBiubkSPzvrHr\niJXt+h4Y8OLb2o7i24SNYkEAJKDz0ax6gyQHcN/sSVfIqp48eQygBAM29v5xoWNIOHmFuLfWaWe7\nmXnbZvI638KT7Ge2GXQT3IS1BYAQ4nZgEfANKeUNKb+fDPwjrZzAY8ClwK+BfdEB8GEp5RO6NKii\njKT37qOlDWAXs+M6BmnyvQwIJ69QlGqzidD5Fp7Eyrc91pvgJqwlAIQQFwMzpJRvF0JsFEIskFLu\nSdz2fmC9lPI7QohbaOX//SnwFSnlNdXILocykr4JH80Kc87LC+wKqu/lyQzSGJq6ooFGf4tax7on\n31zXCHwusCk6HwTOTt4gpfyClPI70b9zgAPAUuBCIcQjQojbhRDTBJAQ4kohxJAQYujgwYOa5E3F\nyv55XHvBG/yZ1VeFDQOZD6kJVd+rAwyEU+BD2+vC1LdooNG/kpHXk2+uKwD6gOei81FgbtaNQoi3\nA6dIKXcCjwLLpZRLgOOAaWEJpZQbpJT9Usr+OXPmaJI3FWn5dJtgoXeKhuQwNYoqTMckw7LR9i4Z\nqom6PGGIZaDi/ZfJZzwZb7oC4BBwfHR+QlY5QohZwOeBP4kuPS6lfCE6HwIWaNZfGXkfr5HCIcvF\nUBWdNqtWQRWmY5JhqbR9WSbrkqGaqEuVIRoUbFrjPFa/imYhk894Mt50jcCP0VL77ATOAoaTNwgh\neoDNwLVSyh9Hl78khFgL/AC4CFinWX9l5BmA2h9t574Rbl61uBmJ3m3FcO9kVNFfu9Z9lzW2uqTP\nRF2qNqeqQehi0PICitU/651XFz7nnVNJEiq7xZIHcBLwfWA98ENaQuCGxD1XAS8D26PjUuBNwOPA\nE8DaonrqzAdw+caHmxnq13QM9wA/4MPO2Lr72KGXpNy6Vsqt66rns5aaO/hTvoOPoSWwHQoCOAVY\nBbxGt4yiw4UAyPp4jQ3vYBtNorUKOv09q+RxMBycrXT9NkJmVEDdweXSoCoAtPcBSClfZtITqLHI\nWgbadhEbHZvgngcfZ+XM79LLr2DHTa0fTKptbLiaNdUHvyxsvKcnrn+A3vvVsZEtrX7VkBmOkKbm\nUd17ZHs3chG6PhaQCRdRHWPS5qH9PLf97+ndcT0gsg1gVYxeusa5vDqLjHUNdOdLxcIVsOA9xflz\nyyDre9TRZjpeKG3GOmdBPQZMTwynSaR5GarGB0u7z6UTSleEgsiTsiZm+jrGpJX987hn4gOMz1xI\n75LLsjt1lZmornEur86i2VVTVwjJ2fnT32wZ008729yu3azEPXW0WR0b/VRWQD6tktCfoasaf5P3\nOc8WpqInquswZQMwoaPLswlYtRfUoYuuUmdTdeeWE8hbr6tK8hVXaGBuCdf6/XZ9l298uBI/IeQD\nmIQJBn3r9r1y8TVflg/d9an6B5ItmGIUvjCcMnBJs4268hinL0zVl+xyJeDaw8dUfUEAGMbIocMt\n5l91IHnWwafAFKPwheH4CFvf3+QKwOc+agIV389Ht88kVAVAV9gATGBWXw9LL74ads2u5v3gs448\naTPQ1cc2OECYMnTbxtb3z9Ppl9X3+9xHTaDi+zUhzr8qggCIo2hQ981mdPFVkVHoRD23rTTm6Ivh\nK8kodAeKJ+F+rUK3bXSFo8s+0ukCvOL7eb+7twS63g10ChTcJiunf0xzZfM1EFbcVbBT3DuTiL9X\nmXfUDeal68roso946m5pDAXvV+SGmeb2qQMfYo519AqgtAuXwszAivT3ZcaVlxPAkzR+SigzW47P\n5EH9HU2tclRp9aWP+IKxEXhkAyBhyRqjwsqViscHVVJHC4AyDdwSFj9nZf9VucLCyg7hMszEpiog\nT61RlQG5VGGUUc+kvVedWdCykNdHfFEhusSuu2HHZ1vnPX1GJyWuVDw+qJI6WgCUaWAfpLESbBro\ndBOIqzAgl4bFMsIq+V6+ZkHLQ1rbdrpQWLwaJsYBaVxgm57kZWkifMg+2NECoEwD1y6NfVAF6Ko1\nVJi7Dt26TKxJRmhdWuNtk9a2ne7J0zcbll1bNxWFcL6ztyQ6WgCUQe3SWFcV4MNMT4W56zC6Tmdi\nRRgbgUduAwQsuTLdcQBabWNaZddU+DAeYtg8tJ9twwdZduYcL72GggAwDO3ofhXj9uzcN8LCi6+r\nJ3mNrRl3tzKxNnbdPRkltqd3ahsXtU2Nq6BaI1x6NmmIaxZU2sJ12wUBYBjatgTdAbt4NTv3jXDV\nk2/kT+fv926JWQlNUuXoIm/Gung1TIwBYjqj97htarWneTZpSGoW2gx++aK53Lf7xWmM3nXbBQFg\nCPEPC5ZtCQmmsfDi6/jT+fu9XGIGFKAo8uqyT7inqSJqtad5LBhharrZNLuA67bTFgBCiNuBRcA3\npJQ3qN6j8lwTYU1yp80QE0yjdvtFwHT4YNSvCVb648E9MPgJeM+6Vj6ChqLN2JcvmsvS01+cxuhd\nj2UtASCEuBiYIaV8uxBioxBigZRyT9E9wJuLnmsqrEnutBli05mGZ4Y6KzDh3990mPzOg59o5WcA\neP/m6rTVhDiDP+OcE2qmRj8UxLlMpoMcBM5WvKfwOSHElUKIISHE0MGDBzXJcw9T28OnIS3kQF8U\nkG7X3dPDFvgWsiGNHl9DX5iEbqgIF9+vah2qz5v8zu9Z18rQ9p511cuqCBshHOoKC6ErAPqA56Lz\nUWCu4j2Fz0kpN0gp+6WU/XPmzNEkr4OQFbcka3DVyVxVmb0uc8yqwzCMDMakkHbFNFXqqVpH3vPx\n+qt85yTmLGjN/D1Q/1SOBxZDu6/d+eCzxsosA10bwCHg+Oj8BNIFSdo9Ks8FqCBLDVSnekhVXVVF\n7eHAzc+YPUcnzlDV72drU57q80X7E8qipBop043SoDrKpLq33deuPu+3Kucm14GuAHiMlvpmJ3AW\nMKx4z08VnvMeyU5Wq99zEnXqlE0ze9U6DNsUjA1wnThDVXcGt5PYR/Wk9s2q3yTv+cQ7Vx4bJQV+\npvBWLEeFXpOG2pX98xifOAKIWviHrgD4GnC/EOK1wAXAe4UQN0gpr8u5ZykgU641DslOljdjtCYc\nPNvwArgRPml1GG4LYwPcZZyhjDZw7pOfeOfK9ZdcrWQKb8VyXLfXrL4eentmcuO3nqK3Z4Z7bz6V\ntGFpB3AKsAp4TZl7VJ5rHz6lhIwjnhJu5NBhuX7wKbl+cDg1RZy1pNK6ae06Md1f8p2y3tHTdzeS\nYjDj3epOX1h3/WVRB7026iTkBHaDIgbv3QDohny9We/o6btbmySUgafCMUAPqgIg7ASuiCJ9sfGN\nHVX13U3fQ6ACHw3kOahkczBl/9BVo3XYng6b9jyvbIURghdORVjz/89CVRe+Kun+xkZg242wbV2h\nO2Ot6e6y3rHOVIc5qScr9SHGKNcUAAAf8ElEQVRTbr+6LpsdtqfDpIuny7J1EVYAmnAuzTO8PJxC\nNQvT2AhPf3WAW598I+BfDPRaEJ9hT4y32nFivHpMe1OrmgrBCCvVb3MFkVV2Tp02Y/HUnnMkBUEA\nqCClw6h6CxgTFD54/eRlYYq30a67WfrMAF977dmctOjcOij1D3FG+cD/aZ3/+1j1cusOJVG1fpv9\nOqvsnDpNqGx9zgCWRBAAKkjpMKrS3JhbmQ/667wsTPE2WrwantnGqfu2wQ/uaETmJsDubDTOKI/r\nm/q3DHzVuevSldWv4yvep7+p97412YIak16WIADUkNJhVKV5rqAYG4FHNgASlqzJ7+CuZnomBnLf\nbJi3BPZto7X1oyFwtcpacmUrwYsOA/JhJZgGXbqy+nW7vGe/NxkELu2+vP6aVbblseSjqicTKq5C\ndR2m3EB1XTFbPv7Dcv3gU3bcONtuiT65Jqq6Sha5DTbMrXDk0GF5x+CjcmzrzfXTnNd2ab9VbWud\n51X3XuiiXd6BpyfLTavDU9deHZh0GSe4gU5CZUmWprfbPLSfgS2taNW9PTPNL+fydOp1oWhJ3p5p\nFc346tZNl8Tmof3cuOVFfnXBf2NN3+x6VS1FSWJM74TWeT75TJwuE20XL29O9PeBAfOh0T1SqdWh\nOuoKAaCyJEtr/FacjqOAtLOcy9Opu0JyABQtyaH1u8bA88IPOm3Aj43wR7/+Z/7jee/m99vfuU5V\nS9m2tRncTecZW20XrzOeIMZXI3RJ1KI6Ulkm1HW43Ans3Y5dVzCl8lGA8R2vOjSlvW/aNRMqjYap\nwY6hDpVSWdx9Seub3X1JtXIM0Ooj7yCogMrBRxctJ0jO5LKWxAZUOsZnODqzt7wInabVcB7NLkuh\nKt0u1H/txDBVE8QYoDVLdePFircIKlKirqMJsYCK4OPsIBft2fDWdfozI1czX8P1TPlWKisjjw3h\nlfpdU1cuNSGrreuM8US3rwDqlr7t+scnjjKwZQ87941w86rF/s4E2mjPgifG9GeBRTNIU4a35Oyt\nYrlTZnL9CquCtPdM0lDTzL+SQbFhBvy6kaU9aII7aMfGAmqnWLvzwWdrqX9yAEqWnTmHbcMHvYoB\nkon24F+yZjI2TNk0jAtXtPK3tsNWJGErfoxKuTnvsrJ/HtefN5c/+vU/ty4UxQ1Ki59T9t0spbhc\n2T+Pay94A6sW9fqVI7qLkBfjKTdWlsO83h27ApjcgDR1I5KrlUFS+rfrbAzis8A097s8PP3N1uad\n086edOGLw5bOPa/c9sx8Ygx23NS6lhIG4LLeB+A710PPjNbvZTcalX03S3aCY7PStG/nketjtyJ3\nhebQdlRaAAghbgcWAd+QUt6Qcc/JwD8CM4Ax4FLg18C+6AD4sJTyCR2iVXDZO15Pb8/MaUzXla9t\nclnYaAOzabdEWyqGvHLbg+qcj7dm7QtXtJhjkgkmaS87GMu+m+0QH2nlN9U4HUfDhVg7FeT4xFFG\nxyamTkZdhn1RMRS0D+Bi4I7ofCOwIOO+DwLnR+e3AH8AvAW4qUx9NozAqsaxvPsaZ9j1AXUaFg+9\n1DJqb107Wb9D91fv0Anv1AE7gG0aibFkBD4X2BSdD9JK8L4nRah8IfbvHOAArfy/FwohlgFPAGuk\nlEdK1l8Zqu6eeSuFJgV78gYqhmHVuEg6de/4bGvm3y5XdZaVNaO3Sa9tdIKR14fgiBXhg5E4VwAI\nIW4DzoxdOge4PTofpTWrz3v+7cApUsqdQoijwHIp5QtCiLuAFcDXU565ErgSYP78+arvYRx5Hyf5\nW90eR41A0YBVzTVgqm4TYYxt0RtQDFtCzKFqyYe9R7kCQEq5Jv6/EGIAOD769wRyvIiEELOAzwN/\nGF16XEp5ODofAhZk1LkB2ADQ399fSyjJIoae/HDOcwM0EUUumzbjItlgFibobbgeuyPRCfaREiir\nAnqMltpnJ3AWMJx2kxCiB9gMXCul/HF0+UtCiLXAD4CLgIpb+OyhrIrHeW6ATkBaMLGsuEgqjNI1\nMzURx6kMswnCohgm2qgDVEtlUFYAfA24XwjxWuACYKkQYhHwPinldbH7/gct9dAnhRCfpGUI/gzw\nZUAAX5dS3leZeksoq5szkhug21AmEYgKo9SZudXNVMswmy6bmWp9GxNt1An2kTJQsRTHD+AUYBXw\nmrLPlj1shYJw6cVjra4mBOzSQdu7ox3sq01jEa2mAsO1y0p6DVVBNweW04WOl0+3tVEOsBUKQkr5\nMpOeQI2ES1WMTl1KtoI6YsCroOqsuj0bXriitZFMMaTCKCey+ciFrOREZpWtKzkDN23gDTPT8tBR\nxZRso662yUXo4J3Ak0h+aJeqmOWL5rJz3wjLF81Ve2BshKe/OsCtT74RyBEadcSAV4Eus4sLjmQi\nEAVoCfUsw/TCFeoGXhWBl2jrjmA8ttVnDgResMl1iQBIfmiX7lf37X6RbcMHWXr6i5xxzgnFD+y6\nm6XPDHDLb1/Nwv4/zL6v6gBRmVXrMKoUweJiRWNEqMdpUDXwqtCdaOuOYDwdYJMINrkuEQDG/fZL\nzH5Kd7KIcS5dvBpqnh0amVWrllNxRWJEqJvOjJWBvD7RmNWBL94yFVYiPvjh1w4VQ0Fdhy0jcOUt\n2E3Yhm4g0bcpA3YInaEOrb7ZzUZmQ2Ox0/oo3Z4PALJnU5WXfhk63eWL5nLf7hf9mL0ZSPRtaoZU\nupy63TNdIOMdtfqmCXWMjyodDfuKCtL4Qkeo5TTQ0QIg66NWZmwZOt2d+0bYNnxwWn21wJE6wwp8\nZEamkfGOWn3TxHczUYZpwa1hX1FBGl/oVntARwsAVx+1Hdr13yZ+ze/8p9+ovxPFByKkhz1Ogy+u\nhr4IIpsw+Y6J76ZlR1D59lkMXiHXghYs9YM0vtCt9oCOzQgG+Rl5TNfT2zOTDffvo7dnhr36kpmC\nsjIHxbNS2cq+lUdXVbSZUd9sp9mRnCL+jobRnuEaz0CX1ZeOzdTF9AxpVWCpjeLMPjMrlwZys3yV\nuL9sOVXQ0SsAl3Cy2kguibOWyGkzJ9VB6WoLvmo93aAOMgxrfTFrRh6/3iCbTVm9f56tb3Rsgo9t\n2jWpAu7/jcL+nVW/U3uEiqW4rqOKF1CjrPqqHhjJ++L/m/LicLUFv5sTshhEYT8P7ZeJsjyi7aH1\n/i8+JE+95l65fvCpab9dvvHhVnkK/TurfhO8C0UvoNqZfN5RRQBkudN5KRgKOosSzaZcU23F3bFR\nRtNg4Z2nMZ4kTPQLz7+VqzHdrmftvbsjATCcTUNOm7mgV1UAdKwKKGsZ7KW7V4Gxy8VGqmNQMQaa\nim3TjvTZMNWBNiyos1b2z+MHe/bxm3s3cs+DH+Cy8/un3mCiX3iuhnOd53t0bILZJ/TkG5H7ZjO6\n+KrIGH/iFLugVzxIRUrUddSZE9jW8zrle7dqqVPdVBWuZrNp9ZisO1bW2Nabpfz0Sa2/NhBWAFrI\nWp35tAKoncnnHbZ2AleBzUTOLsqXUvozoOugw5XQKVtP2baIl++6HWvuP42YJMkWTZdvfLgWVbSq\nAOhYFZAtmPSwSPPXrsWbSAU2dufWse/A1R6DsvWU/SZJzxuX7VizSihNheJCrVJ2f8Wsvh5uXrX4\n2DNxtOkdnzhKb8+M+qIHqEiJug4fVwAqUJXuVWb7lWYQOjO4rWtbM86ta8uX68uKw2e4Mr4bwOiB\n5+VDd31Kjh543mm97T6/98ArtawATK7O2/SuH3zKyoofWysAIcTtwCLgG1LKGzLumQnsiw6AD0sp\nnxBCXA+sAB6RUn5IT2T5D9XZSJXZfmEdeTN2rRmjSPxNQdbM0HMjohdwZXxXQV7fGRth+Nu3cNWT\nb+TPX/Mcl/VucmbEz+vzLnbymlydxw3KvT0z64seoCIl2gdwMXBHdL4RWJBx31uAmxLX3gpsocVB\nPg0sL6qv01cAVuswreuuMkP1ZOaaiybTaJr2vL4T/fbQXZ86ZnyWd1/ipN181PP7CmwYgYHPASui\n8/cCV2Tc90HgSeAR4HZaO44/Cnww+n0pcH3Gs1cCQ8DQ/PnzKzdEGYNRXgdrXOdrAkPzCS49klwy\nbB3k0ZfcfNjO3exICNhG48Z5BlQFQK4KSAhxG3Bm7NI5EUMHGI1m+ml4NJrhvyCEuIuW2qcPeCb2\nbGqORCnlBmADQH9/v8yjTwVlDEZp19uGn/GJIwxs2TvtGW9RZ2C3JoZzdhmAzrQqp4j2st8jr+8k\nf7voVvjan8KewVYdDVfxeeWj7wC5AkBKuSb+vxBiADg++vcEsoPJPS6lPBydDwELgEOKzxpFmt4u\nS5eXdr3dIa4+bwHXXvCGKb9leQV4ndXJBXNOMrgmCARXAnNspJVv+JxrzAdNy4JN20Hf7JYQiEef\nbTBse/n5hrJG4MeAs4GdwFnAcMZ9XxJCrAV+AFwErAMmgFXAP0bPPqtBb2mkGYeyDEZp1+MdIvkR\nk7OFRqwWbDGDOJNPzkhdGoF9Fza77oYdn21FzaxKn+q72l7d+BJG3ADiPKAqA2/EakJFTyQn9fMn\nAd8H1gM/BE6m5RF0Q+K+NwGPA08Aa6NrrwIeAAZoCY7XF9XnuxE4qS9su4mtHxw2pkc0qZMcOXRY\n3jH4aMt4Z1pfm6eHtmiPmNY+vgeZM1lvGd1/sAmVRprbZ5nxWKc9AVs7gYFTaM3kX6Px7PHAJcDp\nKvf7LgCSsPHBTfoep5ZlijHUxGCmvZMqHU3I61yEMm3eCe/rGGnj2clOfQOwJgBcHk0TADZQVqiU\n9mSyzRgsCwZtodvJM2LbcYi6GE3xElIVAB2dEUwVJjLwlCmjzL1ls5rlZYOa1deKYLh5aP9k3YtX\nm83ilIRuRjLFTGDaWd8sZuSqHWlt7vJ9OzWLG+n9zWUGL9MIAgAzKfTKlKF6r07HWtk/b5q3Um7d\nZRlD2cGtK2BcpLIsgkoKzjLt4Yox2hbqRXjktta3e+S2eupvw1F7W0vB6QAhGBzTXb90rP9l3MdU\n79XxIijaEl/Zza2sR4+uh4gPieFVUnCqtsfYyKS/fNG9ZZDmCVS7V45C2BAXcOR95iSAoy2o6Inq\nOuqyAeR585TJ/FMVXuobi963k3TNeSk4s+7JQtvWYnrHrI/GXY0+YKWvZ9Dh5bgyDIINQB9tNQrI\naUu7acs9i6qKPP12lnrIuj6ySGXkg+rGFJLvmvbuKiq0+Oavi27V18OnqTTi6h5fdO8a9gYrapQM\nOnTrMjK2FL+RK7tCUAGlIC9S37TlnmNVRVs99euxEV5+YCP3TExNA1j75hMT7eH7Zq6yMLX5K02l\nEVf3PDDQ2KirLtUounUZGVuKailX4zgIgBwo7SJ2rG9td4y7znyIdx/3FcZnLgQmBUDt+kgT7dFp\n4aNNTRKKyvHBbqKJsuGcq+zSzaqrqMyV/fMYnzjC+MRRRscm9MI7KH4jZ+NYRU9U1+HSBtAUvWCb\nztEDzzdT124r8Ukn2R58g4G2NT2+bGzIUimz0zaChRVABNUll6kAT7rlTJm9zGng7Fhldq+ziui0\nVYNPUGjbov6sq9LIKtfGDFmlTN16fQ0MFwRABJuumTbLcQkjndiWmqJKuQUZsDrKHqGDgrYdHZvg\nY5t2sW34IJDen6vq3ccnjhyzx83q67GSASxu+7ttxzMsXzSX+3a/OKW/xzdTlhkHvo73IABKot2B\nly+ay207ninsBC5nMCaQx+Szop+WEgi2bCZVys2b4fqQhrFuFLTt5qH9bBs+yLIz52T2Z12G3S5v\nfOKoMwba7uc7942kCjUdZu7reA8CIILqR2135Nt2PKN0f1a5LnKYZiKH2eS1Q7IT1zKrscEo82a4\nrgyrDVZh5YVM10V8cjHpkTcjl4FWXaG2n1++qJWravmiuSw9/UWlvCFFqHW850HFUFDX4bMRWPX+\n5H1eGJtzNg85DXerY1z0ceNTWRx6Scqt66TcujZ/k1kXQ8fYWtVA2xQDrwoI0UDdoCwT9KKT+cJs\ndJi5L7RXQfu9my7ILEJngqUzIan6vA7NLqAqAIIKqCLKqkG80AXWHismgo56RYF2I8Zqk6qmZFmL\nV7d2BiMb6bPvAqoqk+T4K6tmqfp8XllNQNcLgCJmobI5JP63CN7qAutAGUFUgiG73LGpVVbfbFh2\nbbUyuxjxMVl1QrV80Vx27hs5pvevAlVafHIJLR0LSAhxuxDiISHEdTn3XCWE2B4du4QQtwkhZgoh\nfhK7/uZqpJtBUVyQot+149ErwHWeAq9RIsZQUUhsJeSEVC7dpmllVY3boxCqumO+fQLxMVl1/N23\n+0W2DR/kvt0vKj+T1a6qtHgVPlpFT9Q+gIuBO6LzjcAChWc+TytWwVuAm8rU58IGUGSk3XvgFXn5\nxofl3gOvpP5etvwyMGEvSJZhWk/pTO/pke7fiB2nqjE7+XxKebbtTXXpvE3nyS5bVtV2bde598Ar\n1toPG0Zg4HPAiuj8vcAVBfe/DtgUnX8QeBJ4BLgdmJnxzJXAEDA0f/786i1RknEkP27R/0UdqEpn\n2XvgFflnXxyUL337r7QZX1bielNMwQujdgG8FHpVBZpCqGrbDLoJ394GTBmObbafEQEA3AZsjx0S\nOCv67T3AxwueXwcsi87/C/Cb0fldwB8UEWdkBVByplW0IijLUKt0llu375VrP7HGqLdIlQ6b9qwX\nng8FzLRbGZUudF2cXcGLPhehSt+y+R62VgADwFI5qQ76RM69rwIeAkT0/3+I/fY/gY8V1WdEABiY\nGZVOtJ6Bsp1l5NBhecfgo3Js682do/qwgQIhrzvQXD9nCiOHDsv1g8Ny/eBTXs1MXcycXbe99r4Z\nyypNVQFQ1gvoMeBsYCdwFjCcc++7gIcjYgC+JIRYC/wAuChaHdhHiqdJlpdIlnU+z6ukjFePjsdQ\nK9Z/f+G9LuCFC2saCtxJdT2vdL2J6nYH3Dy0n4EtewDo7ZlZmgYTKVKz6Gq3C6DdRnn90HXbl+lb\nU2ibea8fO79VpIScnLmfBHwfWA/8EDgZWATckHLvOuDi2P9vAh4HngDWqtRnywicJbWzZhYmN4vU\nhSy6i66rGKqa2iZF6IQVQFVD48ihw/LyjQ8fGxdF76a6WrbVRnW3fR58XAGUEgCtcjkFWAW8puyz\nZQ/XO4FVOo+3apACZNFddD0++MuWHVA/TIVHuHzjw8fGR155oS/4AVUBUHojmJTyZWBTpWWHY6gu\nYVWWc76oQcouy7PoTrs+OjbB+MQRrj5vAX+w+LWpAbFUyjZJf9XnuhVV+2sy0FtReb6MjwBFqEiJ\nug5TK4BOnJWU9u8vseR00V66dXTSt/RZXRFgEIdekmNbb5Z3DD7q7FsTYgFNwudZie6MNpmftND4\nlRbaICO8gosZve438flbxqHSHq4Mlk1cNTWR5kzsupveHdfz3L//dzb3/C+vQsGUDgXRROhsF7ex\njT6tTN1t4bP6eujtmcnAlj3HBkpe+IOXF65i5xlX8/LCVZMXM8Ir6LRX2ffQ2sI/NsKsXbewalEv\nm4f2T/s2PoU+UGkPIyErKH7vIlp8arc2vAqXUBWLVzN+zqd53bkf8G7i0hUrAB3YmJ2llVllRpvU\nz+bRuWn3ODc++TaunT/OmnOiiwaTnTiZmUcCa3jfCDc++TZAM1OTQmC5qjPQovYwOcMteu8iWpy7\nrSq0f1NWekrom03vso9yGZPC1puVjYqeqK6jznwANvSzde6w9E3frEVPZMcYPfB8KZfWaVDYHW7S\n1pBGl+3yXT5fGoaS+vjWp1Wg4mpuAoSEMP5B9SN3kqEzC7W+o4JB3OSATHtXU7vLbaCMr38VQV7V\nB76J46TsHiRdqAqAoAJyiHsefJyR7X/PPRMfiHb4pqO97F21qLcV4tfHROEVUesSP2V3eFIlYzJv\nQ9q75pXvw07ivPor7+hNtL9KTo47H/wRILjsHacdu6eJaqKs795+h+WL5jpVEQUBoAgTOtuVM79L\n73FfYXzmQvLCOxzrJA8MONsu7trrwofEOPF3vvPBHzGwZS/jE0f4yPlnGq2n7LvWzdh0fP2r0Koi\ncAa27AWgt2fGsXt86EM6SBtr7Xe5bcczToV/EACKMDEr611yGfTMoFfV6GrQSFuEumedJlBWiE2d\nyQoAHvvxzxkdm3BmoMtjBrbryUJR/cnf15xzRiXjZpHAWb5oLt99+iC//bqTGzXbz0LeWHMt/IMA\nUISRD1M2F6/D3L11zzpNoGp+5sd/+nO2DR9k89B+Z0LQleC1XU+V8osEzn27X+SBZ0Z498I5fnjO\nVETeWHO+qlExFNR1dJoR2ATqNhBWQdUwxSrlm/SGqRL4zBXNpuuxGQCvqcH1XMPE+6JoBO6KjWCd\nhCZvkGmHKR7YsreQfp3NSVXzwyafT7Z1nKbRsQk+tmlX5W9hM6e0Tj2bh/Zz67ce4emv3lAqX7FK\n+VU2PbpoI1/gcowHFVDD0GRVTSt8xVFAFtLvg00i2dZJ75dtwwdZduacSt/Ct5AHK/vncdZP7mTp\nMwOwy6wKssl9Nw9p37DKd3XZTkEANAxN9XyAFu0fOX+h0r0uB0HWYE22dZKm8YkjtI3Hupgm6BR2\nydrErL4ell58dYv5azofqLZn3r1NQtpkxaZNxCSCCkgTruKnmKjHx1gvRXC57FddcsdpSsZi0sXK\n/nlcfd4CxieOtL5PFO5iaPNnGd+2vpQapghF/eDY75zYmvkrCKCq8a3qVGm+fPAFdn7pL3n54AuV\nykmL6VQU58mXMam1AhBCzAX+SUr5rpx7jgO+CswCbpdSbky7plO/DSjPRKIZ2j3j7+TGLS8CZlUU\nSTpMqEJ8UKf4DJXVRlr/MLFKaZUlj/m5f+Qdq9m5b4SHntpP/7Mb2PnTn7Pw4uuMCMIyG7y00hxq\nxLeqUy00/O1bWPrMADu/DUv/6DPa5aTN2Itm8d6MSRVLcfyglRHs28C/Ftz3UeB/R+ffBE5Mu5ZX\nRhUvoLKWdOWt2FEck7GtN1vxTCgd518BPnlR+ERLGSQzY5nE+sFheeo198r1g8NSylYb3TH4qNyx\n8ZNy8TVfNhYewIZXU93fs0r9oweelw/d9Sk5euB5C5Tlw3YoECymhDyJVi7g7QX3fR1YFJ1/HFiW\ndi2vjCoCoGxsDeVGt5zLs+4BZRtNjN8i5fTcuKbLrhTcrmI9TUZT+1MeTLyTqgAoVAEJIW4D4nvj\nt0opPyNEofGrD3guOh8F5mZcS9Z3JXAlwPz584vqyETZpaWy4cXy5qwmG3lV0FRPkFl9Pdy8avEx\nNZDpstO+uem+4I3awSCa2p/y4PKdREtYaDwoxHYp5bk5v/8zsEZK+TMhxEeBnwGXJq9JKb+cVUZ/\nf78cGhrSoi8gIGAqOsHjJkANQojHpJTZAcci2PQCegw4Ozo/C3g241pAQIADdNuGqoBiGNkHIIT4\nXVq6/b+NXb4T+KYQ4l3AIuBhWuqf5LWAgICAgBqgvQKIq3+klFsTzB8p5Y+B84EHgOVSyqNp13Tr\nDwgICAioBqs7gaWUzwObiq4FBAQEBLhH2AkcEBAQ0KUIAiAgICCgSxEEQEBAQECXIgiAgICAgC6F\n9kYwFxBCHAR+XKGIVwMvGSLHJAJd5RDoKodAVzl0Il2nSinnFN3ktQCoCiHEkMpuONcIdJVDoKsc\nAl3l0M10BRVQQEBAQJciCICAgICALkWnC4ANdROQgUBXOQS6yiHQVQ5dS1dH2wACAgICArLR6SuA\ngICAgIAMBAEQEBAQ0KVovAAQQswVQtxfcM9xQoh7hBAPCCH+JOuaQZpuF0I8JIS4Lueeq4QQ26Nj\nlxDiNiHETCHET2LX31wDXak0CCGuF0I8KoT4O5M0laTtZCHEt4QQg0KI/yeE6LHVZor0TLtH5Tmb\ndLlso5J0edmv6hiHsbpzeZcLvtVoASCEOIVW3oG+gls/DDwmpXwncIkQ4sSMayZouhiYIaV8O3C6\nEGJB2n1SyluklOdGYbXvB74I/A7wlfZ1KeUTJmgqQ1caDUKIt9JK5LMEOCCEWG6KrpK0vR9YL6V8\nD60Mc7+XRq8LetLuKfEe1ujCURtp0OVlv3I9DmO0qfAu63yr0QIAOEorzeQvC+47l8kQ1N8F+jOu\nmUC83EEmM6ClQgjxOmCulHIIWApcKIR4JJq5mAzXrUpXGg3nAP83Sjb9L8C7DNKlTJuU8gtSyu9E\n/84BDmTQ64KetHtUnrNKl8M2KkVXBg1e9CtwOg7bUOFd52KZbzVKAETLs/aybDvw51LKXyg8qpWg\nXpOmD5cs90PALdH5o7QS5SwBjgNW6NBUka40Goy0lQHa2s+/HThFSrkzg96qUHlfa32qIl2AkzYq\nS5f1fqVJVxtWxmEWpJS/VOBd1vuY1YQwpiGlXKP56CHgeOAXwAnR/2nXKtMkhBiIyiUqN1PICiFe\nBSwDPhldelxKeTg6HwK0VQgV6Eqjod1WRc/apg0hxCzg88Af5tBbFSrvm3aP0XbSpMtVG5Wly3q/\n0qTL6jisCGt8q41GrQAqwGWC+jLlvgt4WE5uxviSEOIsIcQM4CLg+4ZoKkNXGg222qoUbUKIHmAz\ncG2UXjSLXhf0uOxTynQ5bKNSdGXQUHt7RXA5DsvAfh+TUjb+ALbHzn8X+LPE76cCTwIDtJZ3M9Ku\nGaLlJFodZj3wQ+BkYBFwQ8q964CLY/+/CXgceAJYa7iNlOhKo4HWROGBqK2GgdfXRNtVwMvA9ui4\n1EabpdBzVgotaTRPu2a5ndLoctJGGnT50K+m0RXd52wcptS9PfpbC9/qmp3AQojX0pKc/yIj3Vva\nNUN1nQKcD3xXSvkzU+VWRRW6hBDHA/8V+Fcp5T6faLMBFXrS7rH9Hr61Uxu6dHVbvyoL23yrawRA\nQEBAQMBUdIsNICAgICAggSAAAgICAroUQQAEBAQEdCmCAAgICAjoUgQBEBAQENCl+P94Aex+9SBZ\nkQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x28130e95e80>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签\n",
    "plt.rcParams['axes.unicode_minus']=False #用来正常显示负号\n",
    "\n",
    "#产生n组点\n",
    "def generate(n):\n",
    "    data=[]\n",
    "    for i in range(n):\n",
    "        x=np.random.uniform(-1,1)\n",
    "        y=np.random.uniform(-1,1)\n",
    "        flag=np.sign(x*x+y*y-0.6)\n",
    "        p=np.random.random()\n",
    "        if (p<0.1):\n",
    "            flag*=-1\n",
    "        data.append([x,y,flag])\n",
    "    return data\n",
    "        \n",
    "data=generate(1000)\n",
    "\n",
    "x1=[i[0] for i in data if i[-1]>0]\n",
    "y1=[i[1] for i in data if i[-1]>0]\n",
    "x2=[i[0] for i in data if i[-1]<0]\n",
    "y2=[i[1] for i in data if i[-1]<0]\n",
    "\n",
    "plt.scatter(x1,y1,s=1)\n",
    "plt.scatter(x2,y2,s=1)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 13\n",
    "\n",
    "这题直接对数据做回归，需要模拟1000次。这里利用了公式来求回归的结果\n",
    "$$\n",
    "w=(X^TX)^{-1}X^Ty\n",
    "$$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.50336\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEFCAYAAAD36MwKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAFHZJREFUeJzt3X2Q5VV95/H3xwESngSUzvBgyIQN\nq5VVULd5XKh0XMAFTGBdN1hxkzXqTtYy1Kb2oRw1m4ArlmE37rJEsMbAQthSS42luEBQV6aYCBGa\n4Wl8ICoM4gjaiGHEELKhvvvH/Y007e3u27dv9505835V3erfPb+Hc35n7nz69Pnd372pKiRJ7Xje\nuBsgSRotg12SGmOwS1JjDHZJaozBLkmNMdglqTEG+24qyRuTPJXk0VmPn521/pNJTlnB+gc6fpIL\nk1y4Uu0YsP5HknwnydHjasfuJsmJSb6W5HtJ/tO426OlMdh3b5+uqsNmPR7euaKqXltVt65Uxf2O\nn2RdkjeuVJ1JXp7kvCVsfzDw74GjgV8AHhmy3hU9rwXqPTjJ7652vZ0NwMVV9TPAVWNqg4ZksGuU\n1gFvXMHjvxwYONiBg4HHq+qpqvqbqnpqyHrXsbLnNZ+DgXEF+yHAwwBV9f0xtUFDMtgblWRTkqlZ\nzy9M8t+T/HmS7ye5dIF990rygyTPS7Klm/b5nSSXLHD8zcAngVO6aaEPzTrkQXPr7Y79R0m2J7kn\nyfFd+dWzR8dJqvu5DbgUOL87/u8vcv4fBu4Afrbb/s9nrXtnkm8leSjJr8xqz8ZuyuYbSc5c6Lzm\nTjEl2ZZk3azlVyW5Ock1s7b5rSQPdFND/2YZ7d+U5HVJPpXkC7PK/213XtuTvL0rm+q2/6MkjyW5\nJcm+3brfTvJwN93ynq7s3CSPAqcAn+zqfnG37sxueubhJO+cVe/VSd6a5KokX59V9n+6tvzXJA8m\n2bjQOWuEqsrHbvigN4J8Cni0e2yes34TMDXr+YXAE8DxwJHA08ALFjj+7fSmL64DLgH+J/D6+Y7f\nlU0Bm+aU9a0XeAvwf4GfBl4FbAN+CrgaeOOs/WvOOV+9hD5aB2ybU3ZWV+/+wIvpTc/sDZwIfALY\nCzgJuH2A87pw1vNtwLpZy1uAU4EDu7J/BNzXnftaYDuwdqntn9X3fwX8KnBQV/bTwBeBI4D9gBng\ngK7tTwO/053n3cA/7/bZAby02/fjO9s6z+vnhfReZ8cCBwH3AGd3664GvgW8eedrqiv7k67e++n9\ntfXAuP/f7CmPvdDu7NNV9folbP+ZqroDoBuVPR94fJ5tt9ALwTuBY4AJ4I+HbGe/es8CPlRVfwt8\nIckTwMtm75QkQ9a3kNPp/ZL5Zvd8P+CIqvpSkvcD7wbOoHe+SzG3rZdU1V/Mev4qenP9X+me70vv\nF8t3l1jPTldV1XU7n1TV3yb5TeA3gNPo/QI5tFv9XeADVVVJ7qEXzAB/AVwMfAp4a1X9cIH6TgHu\nrqp7oTciB84GbujW31BVV87Z5zbgGXqDhL/GGYJVY0fvWb45a3mxT3/bQm8+++v0RrGHd8ujrHfu\n8tw2HTlkfQsJvYuCh1XVYcBRwPYkbwA+ANxFb5Q5+AGTveiNwmf7yz71/umsel/UZ5uleM6+Sf4B\ncAu9X9T/gW5+vPNgdcNontvHv0rvL7EXA1uTLPbLbKF/r37n8sycn1olBvueZSkf5bmF3ijzr4C/\nB34wKxzm8xjwoiRrkhySZM0C9d4IvDnJTyX5JXoXCrfSmx7Y+bbNC/oc/+cAkhzKcD4P/FqS5yc5\ngt4vnYOBk4GbgE8D5w5wXrPbuZ7eNNJCvgCcleSwJAfSm8r4xUX2+T7wwiT7dY99F9j2FfSmgK6i\nF9QvmrXuJ/o/yX70+nsL8PvAk/Sm3uZzK/DyJC9N8nzgX9P7N9QuyGDfvZ2b576P/V+M8Nj3Af+P\nXrB/nd5IdkFVtZVecG6nFxr7LLD5VV0dD9AbNf7Lqnoa+BDw693Fwm1z9rkJeCLJd7t6lqyqbqA3\n9bCV3pz0BVX1GL054V8DHuzaPZHkoAXO66PAZJLP0hutP7RIvVuB/0JveuIr9KZG7l5knx8Cf0jv\nl8+D9P5qms/O/vgu8Ppu+3+4wLH/Bric3r/Bw8BmelMm823/fXph/gngy8DHu77ULiiLD8IkSbsT\nR+yS1BiDXZIaY7BLUmMMdklqzFhuUDr00ENr3bp146haknZbd95552NVtejNc4sGe/eWr48Ca4Af\nAecDV9B7D+71VbXzMyaunFs2n3Xr1jE9Pb3oSUiSnpVkwbfV7jTIVMwbgPdX1Zn0Pivi9cCaqjoZ\nODrJMUleO7ds2IZLkpZn0RF7VV0+6+kE8K+A/9E9/yy9Dzp6BfCxOWXPuf08yXp6d+hx1FFHLavR\nkqT5DXzxNMnJPPsZzdu74sfp3XW3f5+y56iqjVU1WVWTExNL/XwlSdKgBgr2JC8ALgPeRO8zJXZ+\nZsUB3TH6lUmSxmDRAE6yD73Pan5HVT1E72NcT+1WH0fv8zz6lUmSxmCQtzu+GXgl8K4k7wL+F/Ab\n3SfjnUXvSwkK2DynTJI0BoNcPL2C3tsbfyzJdfS+jOCSqnqiK5uaWyZJWn1D3aBUVT/g2XfBzFsm\nSVp9XuSUpMb4nada0LoN14+l3m3vO2cs9UotcMQuSY1xxK5d0rj+UgD/WtDuzxG7JDXGYJekxhjs\nktQYg12SGmOwS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLUGINdkhpjsEtSYwYK9iRrk2zu\nli9Ksql7fC3JO5IcmeTbs8onVrbZkqT5LPqxvUkOAa4B9geoqj+Yte4TwJ8CJwIXd9+PKkkao0FG\n7M8A5wM7ZhcmOR74dlVtB04C3pJkS5L3jr6ZkqRBLRrsVbWjqp7os+rfAZd1yzcCU8DxwMlJjp27\ncZL1SaaTTM/MzCyjyZKkhQx18TTJwcDPVNU3u6Jbq+qHVfUMcBdwzNx9qmpjVU1W1eTEhFPwkrRS\nhn1XzLnADbOe35Tk8CT7AWcCW5fdMknSUIYN9lcDt8x6fhFwM/CXwAer6v7lNkySNJyBv8y6qqZm\nLf/6nHU3Ay8ZXbMkScPyBiVJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLUGINdkhpjsEtSYwx2SWqM\nwS5JjTHYJakxBrskNcZgl6TGGOyS1BiDXZIaM/DnsUtq07oN14+t7m3vO2dsdbfMEbskNcZgl6TG\nGOyS1BiDXZIaM1CwJ1mbZHO3fGSSbyfZ1D0muvIrk9yW5PdWssGSpIUtGuxJDgGuAfbvik4ELq6q\nqe4xk+S1wJqqOhk4OskxK9dkSdJCBhmxPwOcD+zonp8EvCXJliTv7cqmgI91y58FTp17kCTrk0wn\nmZ6ZmVleqyVJ81o02KtqR1U9MavoRnpBfjxwcpJj6Y3mt3frHwfW9jnOxqqarKrJiYmJZTdcktTf\nMDco3VpVTwMkuQs4BngS2LdbfwBelJWksRkmgG9KcniS/YAzga3AnTw7/XIcsG00zZMkLdUwI/aL\ngJuBvwM+WFX3J3kE2JzkCOAsevPwkqQxGDjYq2qq+3kz8JI563YkmQLOAC6ZMycvSVpFI/sQsKr6\nAc++M0aSNCZe5JSkxhjsktQYg12SGmOwS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLUGINd\nkhpjsEtSYwx2SWqMwS5JjRnZx/ZKWp51G64fdxPUCEfsktQYg12SGjNQsCdZm2Rzt3xUkk1JvpBk\nY3qOTPLtrnxTkomVbbYkaT6LzrEnOQS4Bti/K/pt4K1V9dUkNwIvA34BuLiqrlixlkqSBjLIiP0Z\n4HxgB0BVvauqvtqteyHwGHAS8JYkW5K8d0VaKkkayKLBXlU7quqJueVJzge+XFXfAW4EpoDjgZOT\nHNtn+/VJppNMz8zMLL/lkqS+hrp4muRo4D8Cv9sV3VpVP6yqZ4C7gGPm7lNVG6tqsqomJyacgpek\nlbLkYO/m3D8CvGnWSP6mJIcn2Q84E9g6wjZKkpZgmBuUNgBHAZclAfgD4CLgZuDvgA9W1f0ja6Ek\naUkGDvaqmup+vh14e59NXjKiNkmSlsEblCSpMQa7JDXGYJekxhjsktQYg12SGmOwS1Jj/KINSWMz\nri8X2fa+c8ZS72pxxC5JjTHYJakxBrskNcZgl6TGGOyS1BjfFSPNMa53akij4ohdkhpjsEtSYwx2\nSWqMwS5JjTHYJakxBrskNWagYE+yNsnmbnnvJJ9J8sUkb5qvTJI0HosGe5JDgGuA/buiC4A7q+qf\nAK9LcuA8ZZKkMRhkxP4McD6wo3s+BXysW74FmJynTJI0BosGe1XtqKonZhXtD2zvlh8H1s5T9hxJ\n1ieZTjI9MzOzvFZLkuY1zMXTJ4F9u+UDumP0K3uOqtpYVZNVNTkxMTFMWyVJAxgm2O8ETu2WjwO2\nzVMmSRqDYT4E7BrghiSnAb8IfIneNMzcMknSGAw8Yq+qqe7nQ8AZwBeB06vqmX5lK9BWSdIAhvrY\n3qr6Ds++C2beMknS6vPOU0lqjMEuSY0x2CWpMQa7JDXGYJekxhjsktQYg12SGmOwS1JjDHZJaozB\nLkmNMdglqTEGuyQ1xmCXpMYY7JLUGINdkhpjsEtSYwx2SWqMwS5JjRnqq/GSvBU4v3t6MHAnve88\nfaAru6Cq7lt+8yRJSzXsd55eAVwBkOQy4Frgsap6+wjbJkkawrKmYpIcCawFJoHXJLk9yZVJfuIX\nRpL1SaaTTM/MzCynWknSApY7x/42eiP3O4DTq+oEYG/g7LkbVtXGqpqsqsmJiYllVitJms/QwZ7k\necAvA5uAe6vqkW7VNHDM8psmSRrGckbspwFfqqoCrk1yXJI1wHnAPSNpnSRpyYa6eNp5NXBLt/xu\n4MNAgOuq6vPLbZgkaThDB3tVvXPW8lbg2JG0SJK0LN6gJEmNWc5UjFbJug3Xj7sJknYjjtglqTEG\nuyQ1xmCXpMYY7JLUGINdkhpjsEtSYwx2SWqMwS5JjTHYJakxBrskNcZgl6TGGOyS1BiDXZIaY7BL\nUmMMdklqjMEuSY0x2CWpMUv+BqUkewEPdA+AC4DXAWcDt1fV20bXPEnSUg0zYj8W+EhVTVXVFLAP\ncCpwAvC9JKePsH2SpCUaJthPAl6T5PYkVwL/FPizqirgJuC0fjslWZ9kOsn0zMzM8C2WJC1omGC/\nAzi9qk4A9gb2BbZ36x4H1vbbqao2VtVkVU1OTEwM1VhJ0uKWPMcO3FtVT3fL0zwb7gAH4AVZSRqr\nYUL42iTHJVkDnAfsT2+OHeA4YNuI2iZJGsIwI/Z3Ax8GAlwHvAfYnORS4J91D0nSmCw52KtqK713\nxvxY906Yc4BLq+rBEbVNkjSEYUbsP6GqngI+MYpjSZKWxwudktQYg12SGmOwS1JjDHZJaozBLkmN\nMdglqTEGuyQ1xmCXpMYY7JLUGINdkhpjsEtSYwx2SWqMwS5JjTHYJakxBrskNcZgl6TGGOyS1Jgl\nf4NSkoOAjwJrgB8B5wPfAB7oNrmgqu4bWQslSUsyzIj9DcD7q+pM4FFgA/CRqprqHoa6JI3RkoO9\nqi6vqs91TyeAvwdek+T2JFcmGcn3qEqShjP0HHuSk4FDgM8Bp1fVCcDewNnzbL8+yXSS6ZmZmWGr\nlSQtYqhgT/IC4DLgTcC9VfVIt2oaOKbfPlW1saomq2pyYmJiqMZKkha35GBPsg/wceAdVfUQcG2S\n45KsAc4D7hlxGyVJSzDMiP3NwCuBdyXZBHwZuBa4G7itqj4/uuZJkpZqyRc6q+oK4Io5xReNpjmS\npOXyBiVJaozBLkmNMdglqTEGuyQ1xrtEJe1x1m24fmx1b3vfOStehyN2SWqMI/YlGOdveUkalCN2\nSWqMwS5JjTHYJakxBrskNcZgl6TGGOyS1BiDXZIaY7BLUmN2uxuUvElIkhbmiF2SGmOwS1JjDHZJ\naozBLkmNGWmwJ7kyyW1Jfm+Ux5UkDW5kwZ7ktcCaqjoZODrJMaM6tiRpcKN8u+MU8LFu+bPAqcDX\nd65Msh5Y3z19Msn9I6x7uQ4FHht3I3ZR9k1/9kt/9kt/P+6X/OGyjvNzg2w0ymDfH9jeLT8OvHL2\nyqraCGwcYX0jk2S6qibH3Y5dkX3Tn/3Sn/3S32r3yyjn2J8E9u2WDxjxsSVJAxpl+N5Jb/oF4Dhg\n2wiPLUka0CinYj4FbE5yBHAWcNIIj73Sdskpol2EfdOf/dKf/dLfqvZLqmp0B0sOAc4AbqmqR0d2\nYEnSwEYa7JKk8fMCpyQ1Zo8I9kHviE2yNsldS91vdzVMvyTZK8m3kmzqHi9bndaunsX6Zb4+SHJR\nkjuSfGB1W7w6ltEvd88qO2N1W73ylvD/6PIkv7LU/YbRfLAv8Y7Y/0b3ls3W76Qdtl+AY4GPVNVU\n97hvpdu6mgbsl5/ogyT/mN67wk4Avpfk9FVs9opbRr+8EPjarLLPrWrDV9ig/4+SnAYcVlWfWcp+\nw2o+2Ol/R+xPSPIq4EfAzou+A+23G5tiuH45CXhNktu7Ecdu92Uti5hi8X7p1we/BPxZ9S5a3QSc\nthqNXUVTDNcvJwInJLk1yaeSHLg6zV01UyzSL0n2Bj4EbEty7qD7LceeEOxz74hdO3eDJPsA/xnY\nsJT9dnPD9ssdwOlVdQKwN3D2CrdztQ3y796vD/b41wv9++UB4NVVdQpwL/Bbq9DW1TRIv/wm8BXg\nEnq/5C4YcL+h7QnBPsgdsRuAy6vqr5e43+5s2H65t6oe6ZangaamqBisX/r1ga+X/v3yAPCNOWUt\nGaRfXgFs7N4C/r+BXx5wv6G19uLrZ5A7Yk8H3pZkE/DyJH8y4H67s2H75dokxyVZA5wH3LMKbV1N\ng/RLvz7w9dK/Xy4Gdl4wfB175uvlG8DR3fIk8NCA+w2vqpp+AM+n92J6P/DVrhPfs8D2m+bZ76Bx\nn8su0i8vpfcn9X3AxeM+j3H0S78+oDdI+iJwKXA/8PPjPpddpF8OB74EbKU3z7z3uM9lDP1yIPBx\n4BbgNuDIlc6XPeIGpWHviG39TtrWz29Yy3i97AucA2ypqgdWqn3j4uulv10xX/aIYJekPcmeMMcu\nSXsUg12SGmOwS1JjDHZJaozBLkmN+f9JZaknKQGQ6QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x28130d7bb00>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from numpy.linalg import inv\n",
    "\n",
    "Ein=np.array([])\n",
    "for i in range(1000):\n",
    "    data=generate(1000)\n",
    "    X=np.array([[1]+i[:-1] for i in data])\n",
    "    Y=np.array([i[-1] for i in data])\n",
    "\n",
    "    w=inv(X.T.dot(X)).dot(X.T).dot(Y)\n",
    "\n",
    "    error=np.sum(np.sign(X.dot(w)*Y)<0)/1000\n",
    "    Ein=np.append(Ein,error)\n",
    "\n",
    "print(np.average(Ein))\n",
    "plt.hist(Ein)\n",
    "plt.title('Ein without feature transform')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "所以$E_{in}$的均值约为$0.5$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 14\n",
    "\n",
    "先做特征转换，再重复上题的步骤，画出$\\tilde w_3$的直方图，这里把15题的任务一起做了。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEFCAYAAADt1CyEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAADdVJREFUeJzt3X+sX/Vdx/Hny0KTrlVWwk0DGtIQ\nK4ZEmuEFaSzbNQIZ6zSELWKC7o/NdDHLovEfQfgHZQb3B5Gg4Jo0BjEugy3oFAnFBKROJtyCm9Nl\nm9laXSfhzuJqzYIT3/5xT2253Pae+73n9nvv5/t8JN9wzud7frw/9/T7yuFzzvl+U1VIktr1feMu\nQJK0ugx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfTSAkmuTvJnSZ5KctO465FWKt5HL52S5ALgReDD\nQAGPApdX1WtjLUxagfPGXYC0xlwM3FVVzwAk+VfghwCDXuuWQa+JleSfgXcCvwd8EZgFPlBVP5dk\nA3AzcD7wT+OrUlo5x+g1yV4CLgf+E7gU+JGuDeBXgD8GPlFVb4ynPGkYBr0m2UvATuB7wEbmQ/8l\ngKq6D9gB/EaSHx5bhdIADHpNspeAnwG+3s1fARxLciVAVX2T+Quzl4+nPGkYBr0m2SHgp4CvAkeY\nH755G/CpJG9PcjFwNfDy+EqUVs6LsZpYVfXv3V01XwU2Ay9X1XNJ/gj4MvBfwK9W1bfGWae0Ut5H\nL0mNc+hGkhpn0EtS4wx6SWqcQS9JjRvLXTcXXXRRbd++fRy7lqR169ChQ9+uqqnlrjeWoN++fTuz\ns7Pj2LUkrVtJjoyynkM3ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOH94\nRGvS9tufGNu+D9+7Z2z7llaDZ/SS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0\nktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXFLBn2SC5I8meRAkseTbEzyL0me7V4/1i13d5IXk/z+\n6pctSeqrzxn9bcB9VXUj8ApwO/DJqprpXv+Q5MeB3cA1wKtJrl+9kiVJy7Fk0FfVg1X1dDc7BfwP\n8N4kLyTZn+Q84F3AZ6qqgKeA6xZuJ8neJLNJZufm5gbsgiTpbHqP0SfZBWwFngaur6prgPOB9wCb\ngaPdoseAbQvXr6p9VTVdVdNTU1MrLlyS1E+vnxJMciHwAPA+4JWqer17axbYAZwANnVtW/AirySt\nGX0uxm4EHgPuqKojwCNJdibZANwMfAE4xPwYPcBO4PDqlCtJWq4+Z94fAq4C7kzyLPCPwCPA3wPP\nV9VfAX8DvCPJ/XQXa1enXEnSci05dFNVDwEPLWi+e8Ey/9vdabMHuL+qvjFciZKkleg1Rt9HVX0X\n+PRQ25MkDcOLppLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BL\nUuMMeklqnEEvSY0z6CWpcYN9TbHUiu23PzGW/R6+d89Y9qv2eUYvSY0z6CWpcQa9JDXOoJekxhn0\nktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhq3ZNAnuSDJk0kOJHk8\nycYk+5M8n+Su05Z7S5skafz6nNHfBtxXVTcCrwA/D2yoql3AZUl2JLllYdvqlSxJWo4lf3ikqh48\nbXYK+AXgd7v5A8Bu4B3Aowvavnb6dpLsBfYCXHrppSsqWpLUX+9fmEqyC9gKHAaOds3HgKuAzYu0\nvUlV7QP2AUxPT9fIFeucGtevLUkaTq+LsUkuBB4APgicADZ1b23ptrFYmyRpDehzMXYj8BhwR1Ud\nAQ4xPzQDsJP5M/zF2iRJa0CfoZsPMT8Uc2eSO4E/BH4xySXATcC1QAEHF7RJktaAPhdjHwIeOr0t\nyWeBG4CPV9V3uraZhW2SpPHrfTH2dFX1GqfusjljmyRp/LxoKkmNM+glqXEGvSQ1zqCXpMYZ9JLU\nOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z\n6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUuF5B\nn2RbkoPd9A8m+WaSZ7vXVNe+P8nzSe5azYIlScuzZNAn2Qo8DGzumn4C+FhVzXSvuSS3ABuqahdw\nWZIdq1eyJGk5+pzRvwHcChzv5q8FfinJS0l+u2ubAR7tpg8AuxduJMneJLNJZufm5lZWtSSptyWD\nvqqOV9V3Tmt6kvlgvxrYleRK5s/2j3bvHwO2LbKdfVU1XVXTU1NTKy5cktTPeSOs87dV9TpAkpeB\nHcAJYFP3/ha8yCtJa8YogfxUkouTvA24EfgScIhTwzU7gcPDlCdJWqlRzujvBp4B/hv4g6r6SpJ/\nAw4muQS4iflxfEnSGtA76KtqpvvvM8CPLnjveJIZ4Abg4wvG9CVJYzTKGf2iquo1Tt15I0laI7xo\nKkmNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS\n1DiDXpIaZ9BLUuMMeklqnEEvSY0b7BemJK3M9tufGNu+D9+7Z2z71urzjF6SGmfQS1LjDHpJapxB\nL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWpcr6BPsi3JwW76/CR/nuRzST54pjZJ\n0tqwZNAn2Qo8DGzumj4KHKqqnwTen+T7z9AmSVoD+pzRvwHcChzv5meAR7vp54DpM7S9SZK9SWaT\nzM7Nza2gZEnSciz5NcVVdRwgycmmzcDRbvoYsO0MbQu3sw/YBzA9PV0rKXrSjPPrayWtf6NcjD0B\nbOqmt3TbWKxNkrQGjBLIh4Dd3fRO4PAZ2iRJa8AovzD1MPCXSa4DrgD+jvlhm4VtkqQ1oPcZfVXN\ndP89AtwAfA64vqreWKxtFWqVJI1gpN+MrapvceoumzO2SZLGz4umktQ4g16SGmfQS1LjDHpJapxB\nL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS\n1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mN\nO2+5KyQ5D/h69wL4KPB+4D3AC1X1keHKkySt1Chn9FcCn6yqmaqaATYCu4FrgFeTXD9gfZKkFRol\n6K8F3pvkhST7gZ8GPlNVBTwFXLfYSkn2JplNMjs3Nzd6xZKkZRkl6F8Erq+qa4DzgU3A0e69Y8C2\nxVaqqn1VNV1V01NTUyMVK0lavmWP0QNfrKrXu+lZToU9wBa8wCtJa8ooofxIkp1JNgA3A5uZH6MH\n2AkcHqg2SdIARjmj/03gT4AAnwXuAQ4muR94d/eSJK0Ryw76qvoS83fe/L/uTps9wP1V9Y2BapMk\nDWCUM/q3qKrvAp8eYluSpGENEvSS1rfttz8xlv0evnfPWPY7abxDRpIaZ9BLUuMMeklqnEEvSY0z\n6CWpcQa9JDXO2yuXYVy3oEnSSnhGL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqc\nQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIa59cUSxqbcX319+F794xlv+PiGb0kNc6gl6TGGfSS\n1DiDXpIaZ9BLUuPW5V03/ki3JPU3aNAn2Q9cATxRVfcMuW1JGso4TxbHcWvnYEM3SW4BNlTVLuCy\nJDuG2rYkaXRDntHPAI920weA3cDXTr6ZZC+wt5s9keQrA+57KBcB3x53EefIJPUVJqu/9nUNy++s\naPXLR1lpyKDfDBztpo8BV53+ZlXtA/YNuL/BJZmtqulx13EuTFJfYbL6a1/blWR2lPWGvOvmBLCp\nm94y8LYlSSMaMowPMT9cA7ATODzgtiVJIxpy6OZPgYNJLgFuAq4dcNvnypoeWhrYJPUVJqu/9rVd\nI/U3VTVYBUm2AjcAz1XVK4NtWJI0skGDXpK09njBVJLWmSQXJrkhyUV9lp/YoE+yP8nzSe46yzIX\nJHkyyYEkjyfZeC5rHEqfvnbLbUty8FzVNbSex7TX32I96NnfdX1MT1qqr618VqFXX7cCfwFcAzyT\nZGqpbU5k0C/jKd7bgPuq6kbgFeDd56rGofTta/eP52Hmn4dYd/r0s6Wnt3v2d10f05N6Hrd1/1mF\n3n29Evi1qvoY8BQLnllazEQGPYs/xfsWVfVgVT3dzU4Br65+aYOboUdfgTeAW4Hj56Cm1TDD0v3s\ns8x6McPSfVnvx/SkGZboayOfVejX17+uqs8neSfzZ/XPL7XRdfntlcuV5BO8+dHhdwH7u+m3PMW7\nyPq7gK1V9fnVqXA4o/a1qo53669qfavorE9mL2OZ9WLJvjRwTE/qfdzW02f1DHr1NfMH9VbgNeB7\nS210IoK+qj58+nyS++n5FG+SC4EHgPetWoEDWklf17k+T2a39PR2S31ZSq++rrfP6hn06mvN3y75\nkSS/Bfws8KmzbbTlfxxn0+sp3u6CzmPAHVV15NyUNrhJeWK5Tz9b+lu01JelLNnXRj6r0K+vv57k\nA93s24H/WHKrVTVxL+AHgC8A9wFfBi5g/nv071mw3C8z/79Gz3avW8dd+2r19bTlnx13zQP1c+ci\nx/Mtf4tx172a/V3vx3SZx3bdf1aX0detwNPAc8CDdM9Dne01sQ9MTdJTvJPS1z79bOlv0VJflmJf\nV7jNSQ16SZoUkzpGL0kTw6CXpMYZ9JLUOINekhpn0EtS4/4PI5ooQw6xi/UAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x281313f9908>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ -9.92410543e-01   1.64976493e-03  -6.48205307e-04   1.52070206e-04\n",
      "   1.55620104e+00   1.55888923e+00]\n",
      "w3的平均值0.000152070206015\n"
     ]
    }
   ],
   "source": [
    "W=[]\n",
    "Eout=np.array([])\n",
    "for i in range(1000):\n",
    "    data=generate(1000)\n",
    "    X=np.array([[1]+i[:-1]+[i[0]*i[1],i[0]**2,i[1]**2] for i in data])\n",
    "    Y=np.array([i[-1] for i in data])\n",
    "    \n",
    "    w=inv(X.T.dot(X)).dot(X.T).dot(Y)\n",
    "    \n",
    "    #测试数据\n",
    "    data1=generate(1000)\n",
    "    X1=np.array([[1]+i[:-1]+[i[0]*i[1],i[0]**2,i[1]**2] for i in data1])\n",
    "    Y1=np.array([i[-1] for i in data1])\n",
    "\n",
    "    error=np.sum(np.sign(X1.dot(w)*Y1)<0)/1000\n",
    "    Eout=np.append(Eout,error)\n",
    "    \n",
    "    #记录w\n",
    "    W.append(w)\n",
    "\n",
    "W=np.array(W)\n",
    "w3=np.array([i[3] for i in W])\n",
    "plt.hist(w3)\n",
    "plt.title('w3')\n",
    "plt.show()\n",
    "\n",
    "print(W.mean(axis=0))\n",
    "print(\"w3的平均值\"+str(w3.mean()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 15\n",
    "\n",
    "做出$E_{out}$的直方图"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEFCAYAAAD36MwKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAEjtJREFUeJzt3XuQpFV9xvHvEy66XIQFJouIuIUh\nJkYEzcpFsRwtoAJ4ixKx4iXesl4SUsaYCopl0GgKjVohlhLXQCBWNAV4Q4EARrdYlCi7CIIXRGCR\nrIKLa0BKJAn55Y9+Cc0wl56enunl7PdT1TXdp98+76/PzDx95rz9TqeqkCS141fGXYAkabQMdklq\njMEuSY0x2CWpMQa7JDXGYJekxhjsjUnyqiT3JLmt7/LYraCuKwepI8lZSV41RP+nJPlxkh8l2X+o\nIrdBSQ5N8r0kP0ny5+OuR6NhsLfp81W1d9/l1mE6SXJwkheOoqCqetrUOkbVf5LdgbcA+wO/Bvx4\nyH5WDvOislBJdk/y5qXeb+ck4L1V9avAmWOqQSNmsGs2BwMjCfZF7n93YEtV3VNVv6iqe4bsZyXw\nqhHUM1+7A+MK9uXArQBV9dMx1aARM9i3EUl+JckHk2xKck2Sp3XtD1r6SFLd143AacAJ3XLOO2fp\ne98k3++ub0kymeQDSd7Ut83GJCv7b8/S/75Jvtb19dY5ntcngSuBx3b9/GvffW9P8sMktyR5Xt84\nrOmWbH6Q5OiufR3wGeDpXT8f79pPSXLKdM+ju/6cJF9JcnbfNq9OclO3NPSHC6h/bZLjk3wuyZf7\n2t/QPa9NSf6ia5vstv9gkjuSXJZkWXff65Pc2i23vKdre0GS24CnA5/p9v2E7r6ju+WZW5O8vW+/\nZyV5Y5Izk9zQ1/bFrpa/SXJzkjWzPWctgary0tCF3ozzHuC27rKua38d8G/AI4HnABuBRwBnAa/q\ne3xN6eusAfe7id6M91LgTcD5wGF9928EVk5T61lT2s6iN4N8PPA04M4B9r0S2Dil7Zju+e4MPIHe\n8swOwKHAecD2wGHAN/oeMwmsndLPKcAp0z2P7vpVwBHArl3bbwHXAnsAK7pxWTHf+rv2tcD3gecD\nu3VtjwS+CuwD7ARsBnbpar8X+OPueV4N/G73mLuAJ3WPPff+Wvv2Mdl3e8/u5+bJwG7ANcCxfd+b\nHwKvBfboa/uHbr/X0/sr7KZx/x5s65ftUYs+X1UvndJ2DPDxqvol8OUkdwIH9m+QJAvY57XA84Av\nAr9Jb637miH7OruqbkxyE/CoIfs4kt4Lw43d7Z2Afarq60k+BLwbOAqYmGe/U8fo/VV1ed/t59Bb\n6/9Od3sZvReW2+e5n/udWVXn33+jqn6Z5JXAK4Bn0nsB2au7+3bgI1VVSa6hF8wAlwPvBT4HvLGq\nfj7L/p4OXF1V34LejBw4Friwu//CqjpjymOuAO4DvgH8J64EjJ3fgG1LTbk+9T/APWYBfV9Fb738\nSnqzyXtr+LXuG6H702F4oXdQcO+q2hvYD9iU5GXAR4Bv0ptlDt5hsj29WXi/f59mv//Ut999p9lm\nPh702CSPBy4DtgB/Rrc+3rm5b8z6x+75wN/Re4G5LslcL2az/ZxM91zum/JVY2awbzsuAl6b5BFJ\nnkXvgN119P5Mv/9tiCdOecwdwOMAkuzF7K6iN1v9Pr3wu3aAmmbqfxT/cvRLwEuSPCrJPvReLHYH\nDgcuBj4PvGCaevZNsl2S5Um248Hjs5re8tVsvgwck2TvJLvS+6vliXM85qfAnkl26i7LZtn2KfSW\ngM6kF9T79t33kHFLshO97/NVwDuBu+n9NTWTrwEHJ3lSkkcBf0DvZ0cPIwZ7m16QB7+P/cX0guBa\n4CZ6s7ffq6p7gY8Dv98dtNs4pZ+LgTuT3E4vKGezAfhZVd0B3EAvSOYyn/7npaoupLf0cB29NekT\nu9rOAl4C3AzsCEwk2a17zHVdHZu6x+0I/AuwKskl9F6wbpljv9cBf0VveeI79JZGrp7jMT8H3kfv\nxedm4NGzbH7/ON0OvLTb/tdn6fsXwEfpfe9vBdbRWzKZafuf0gvz84BvA+d2Y6mHkSzsr11J0tbG\nGbskNcZgl6TGGOyS1BiDXZIaM5YTlPbaa69auXLlOHYtSQ9bGzZsuKOq5jypbizBvnLlStavXz+O\nXUvSw1aSWd9uez+XYiSpMQa7JDXGYJekxhjsktQYg12SGmOwS1JjDHZJaozBLkmNMdglqTF+5qk0\nxcqTLhjLfjeeetxY9qv2OGOXpMYY7JLUGINdkhpjsEtSYwx2SWqMwS5JjTHYJakxBrskNcZgl6TG\nGOyS1BiDXZIaY7BLUmMMdklqjMEuSY0x2CWpMQa7JDVmzmBPsluSi5JckuSzSXZM8sMka7vLgd12\n70pyZZKPLH7ZkqSZDDJjfxnwoao6GrgNOAn4VFVNdpdrk/w2cARwCPCTJEcuXsmSpNnMGexV9dGq\nurS7OQH8D/DcJN9IckaS7YFnAZ+uqgIuBp45tZ8kq5OsT7J+8+bNI3wKkqR+A6+xJzkcWA5cChxZ\nVYcAOwDHAjsDm7pNtwArpj6+qtZU1aqqWjUxMbHgwiVJ0xvow6yT7AF8GHgxcFtV3dvdtR44ALgb\nWNa17YIHZSVpbAY5eLojcC7wtqq6BfhEkoOSbAe8ELgG2EBvjR3gIGDj4pQrSZrLIDPr1wJPBU5O\nshb4NvAJ4Grgiqr6EnA58JQkp9EdXF2cciVJc5lzKaaqTgdOn9L8rinb/G/3TpjjgNOq6ubRlShJ\nmo+B1tgHUVX3AOeNqj9J0nA8yClJjTHYJakxBrskNcZgl6TGGOyS1BiDXZIaY7BLUmMMdklqjMEu\nSY0x2CWpMQa7JDXGYJekxhjsktQYg12SGmOwS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLU\nGINdkhpjsEtSYwx2SWqMwS5JjTHYJakxBrskNcZgl6TGGOyS1Jg5gz3JbkkuSnJJks8m2THJGUmu\nSPKOvu0e0iZJWnqDzNhfBnyoqo4GbgNeCmxXVYcD+yc5IMmLprYtXsmSpNlsP9cGVfXRvpsTwMuB\nv+1uXwIcATwFOGdK2w39/SRZDawG2G+//RZUtCRpZgOvsSc5HFgO3Aps6pq3ACuAnadpe5CqWlNV\nq6pq1cTExIKKliTNbM4ZO0CSPYAPAy8G3gIs6+7ahd6Lw93TtElDW3nSBeMuQXrYGuTg6Y7AucDb\nquoWYAO9pRaAg4CNM7RJksZgkBn7a4GnAicnORn4R+AVSfYBjgEOAwpYN6VNkjQGgxw8PR04vb8t\nyfnAUcD7q+rOrm1yapskaekNtMY+VVX9jAfeBTNjmyRp6XmQU5IaY7BLUmMMdklqjMEuSY0x2CWp\nMQa7JDXGYJekxhjsktQYg12SGmOwS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLUGINdkhpj\nsEtSY4b6aDxJo7fypAvGtu+Npx43tn1r9JyxS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLU\nGINdkhpjsEtSYwYK9iQrkqzrrj8myX8kWdtdJrr2M5JckeQdi1mwJGl2cwZ7kuXA2cDOXdOhwHur\narK7bE7yImC7qjoc2D/JAYtXsiRpNoPM2O8DTgDu6m4fBrwuyVVJ/rprmwTO6a5fAhwxtZMkq5Os\nT7J+8+bNC6takjSjOYO9qu6qqjv7mi6iF+RPAw5P8mR6s/lN3f1bgBXT9LOmqlZV1aqJiYkFFy5J\nmt4w/93xa1V1L0CSbwIHAHcDy7r7d8GDspI0NsME8MVJHp1kJ+Bo4DpgAw8svxwEbBxNeZKk+Rpm\nxv4u4CvAfwF/X1XXJ/kxsC7JPsAx9NbhJUljMHCwV9Vk9/UrwG9Mue+uJJPAUcD7p6zJS5KW0Mg+\nQamqfsYD74yRJI2JBzklqTEGuyQ1xmCXpMYY7JLUGINdkhpjsEtSYwx2SWqMwS5JjRnZCUpq08qT\nLhh3CZLmyRm7JDXGYJekxhjsktQYg12SGmOwS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLU\nGINdkhpjsEtSYwx2SWqMwS5JjTHYJakxBrskNcZgl6TGGOyS1BiDXZIaM9CHWSdZAZxXVc9MsgPw\nGWAP4IyqOnO6tkWrWNLIjetDyzeeetxY9tu6OWfsSZYDZwM7d00nAhuq6hnA8Ul2naFNkjQGgyzF\n3AecANzV3Z4EzumuXwasmqHtQZKsTrI+yfrNmzcvoGRJ0mzmDPaququq7uxr2hnY1F3fAqyYoW1q\nP2uqalVVrZqYmFhY1ZKkGQ1z8PRuYFl3fZeuj+naJEljMEwAbwCO6K4fBGycoU2SNAYDvStmirOB\nC5M8E3gi8HV6yzBT2yRJYzDwjL2qJruvtwBHAV8Fjqyq+6ZrW4RaJUkDGGbGTlX9iAfeBTNjmyRp\n6XmQU5IaY7BLUmMMdklqjMEuSY0x2CWpMQa7JDXGYJekxhjsktQYg12SGmOwS1JjDHZJaozBLkmN\nMdglqTEGuyQ1xmCXpMYY7JLUGINdkhpjsEtSYwx2SWqMwS5JjTHYJakxBrskNcZgl6TGGOyS1BiD\nXZIaY7BLUmMMdklqjMEuSY3Zfr4PSLI9cFN3ATgROB44FvhGVf3R6MqTJM3XMDP2JwOfqqrJqpoE\ndgSOAA4BfpLkyBHWJ0map3nP2IHDgOcmeTZwLXA98OmqqiQXA8cAX5r6oCSrgdUA++233/AVb4NW\nnnTBuEuQ9DAyzIz9SuDIqjoE2AFYBmzq7tsCrJjuQVW1pqpWVdWqiYmJoYqVJM1tmBn7t6rq3u76\neh4Id4Bd8ICsJI3VMCH8iSQHJdkOeCGwM701doCDgI0jqk2SNIRhZuzvBj4JBDgfeA+wLslpwO90\nF0nSmMw72KvqOnrvjPl/3TthjgNOq6qbR1SbpMaN640BG089biz7XSrDzNgfoqruAc4bRV+SpIXx\nQKckNcZgl6TGGOyS1BiDXZIaY7BLUmMMdklqjMEuSY0x2CWpMQa7JDXGYJekxhjsktQYg12SGmOw\nS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLUmJF85um2YlwfvCtptMb5u7wUH6TtjF2SGmOw\nS1JjDHZJaozBLkmNMdglqTEGuyQ1xmCXpMYY7JLUGINdkhoz0jNPk5wBPBG4oKreM8q++3kGqCTN\nbGQz9iQvArarqsOB/ZMcMKq+JUmDG+WMfRI4p7t+CXAEcMP9dyZZDazubt6d5PpZ+toLuGOEtbXK\ncRqcYzUYx2kwQ49T3reg/T5ukI1GGew7A5u661uAp/bfWVVrgDWDdJRkfVWtGmFtTXKcBudYDcZx\nGszWPk6jPHh6N7Csu77LiPuWJA1olOG7gd7yC8BBwMYR9i1JGtAol2I+B6xLsg9wDHDYAvoaaMlG\njtM8OFaDcZwGs1WPU6pqdJ0ly4GjgMuq6raRdSxJGthIg12SNH4e4JSkxix5sCc5I8kVSd4x6DZJ\nlie5MMn6JB9bumrHa8CxWpFkXd/tHZJ8IclXk7xmaSodryHHab8ka5N8OcmaJFmaasdnmHHqa39S\nkksXt8KtwwLH6QtJDl7cCue2pME+yNmpM2zzCuCfu/eN7ppkq33/6KgMOFbLgbPpnUNwvxOBDVX1\nDOD4JLsuScFjsoBxej3wxqp6DvBY4MClqHdcFjBOdC96HwJ2WIpax2mB4/Qy4MaqunpJip3FUs/Y\nJ3no2amDbPNT4ElJdqf3S3jrola5dZhk7rG6DzgBuGuGx10GtP4iOMkQ41RVJ1fVd7ube9L+2ZaT\nDPfzBPBq4CuLVtnWZZIhxinJHsAHgZ8lefYi1zinpQ72qWenrhhwm8vpnUr7J8B3u/bWzTlWVXVX\nVd0538c1ZthxAiDJCcC3q+pHi1fiVmGocUqyJ/By4AOLXuHWYdifpz8FzgU+BrwyyfMXtco5LHWw\nD3J26nTb/CXwhqp6N/A9ejOI1g17Ju+2dgbw0M83yf7AW4E3L0JdW5thx+lU4G1V9d+LUtXWZ9hx\negrwke5t3ufQm/mPzVL/0g9ydup02ywHDkyyHXAosC28R3PYM3m3tTOAh3q+3Trpp4DXzDSbb8yw\nPxfPAt6XZC1wcJJF+3fcW4lhx+kHwP7d9VXALaMta56qaskuwKOAa+gdiPkuvYF7zxzb7AYcAnyb\n3qvppcAuS1n3OC6DjFXftmv7rj+uG6vTgCvpHQga+/PZCsfpfcCPgbXd5Vnjfi5b4zgN0t7SZQE/\nT/sAFwJf7TJq13E+jyU/QWmQs1M9g7Vn2HHo/q3DEcDFtQ3MRv15GYzjNJgWxskzTyWpMa0fWJOk\nbY7BLkmNMdglqTEGuyQ1xmCXpMb8H0ojvkd0iAyBAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x28131437898>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.125974\n"
     ]
    }
   ],
   "source": [
    "plt.hist(Eout)\n",
    "plt.title('Eout with feature transform')\n",
    "plt.show()\n",
    "print(Eout.mean())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "所以$E_{out}$的平均值为0.126左右"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 16\n",
    "\n",
    "这题实际上是多元Logistic回归，同课件里的例子，我们要最大化似然函数，这等价于题目中所说的最小化负的对数似然函数，先把似然函数求解出来\n",
    "$$\n",
    "L=\\prod_{j=1}^{N}\\frac{exp(w_{y_j}^T x_j)}{\\sum_{i=1}^K exp(w_i^T x_j)}\n",
    "\\\\ln(L)=\\sum_{j=1}^{N}[ln(exp(w_{y_j}^T x_j))-ln(\\sum_{i=1}^K exp(w_i^T x_j))]=\\sum_{j=1}^{N}[w_{y_j}^T x_j-ln(\\sum_{i=1}^K exp(w_i^T x_j))]\n",
    "$$\n",
    "所以只要最小化\n",
    "$$\n",
    "E_{in}=-\\frac{ln(L)}{N}=\\frac 1 N\\sum_{j=1}^{N}[ln(\\sum_{i=1}^K exp(w_i^T x_j))-w_{y_j}^T x_j]\n",
    "$$\n",
    "即可"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 17\n",
    "\n",
    "求偏导即可\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\frac {\\partial E_{in}}{\\partial w_n}&=\\frac  {\\partial (\\frac 1 N\\sum_{j=1}^{N}[ln(\\sum_{i=1}^K exp(w_i^T x_j))-w_{y_j}^T x_j])}{\\partial w_n}\n",
    "\\\\&=\\frac 1N\\sum_{j=1}^{N}( \\frac{exp(w_n^Tx_j) x_j}{\\sum_{i=1}^K exp(w_i^T x_j)}-[\\![y_j=n]\\!]x_j)\n",
    "\\\\&=\\frac 1N\\sum_{j=1}^{N}(h_n(x_j)-[\\![y_j=n]\\!])x_j\n",
    "\\end{aligned}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 18\n",
    "\n",
    "使用梯度下降公式\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\nabla E_{in}(w)\n",
    "&=-\\frac 1N \\sum_{n=1}^N\\frac{y_nx_n}{1+e^{y_nw^Tx_n}}\n",
    "\\\\&=\\frac 1N \\sum_{n=1}^N-{y_nx_n}\\theta(-y_nw^Tx_n)\n",
    "\\end{aligned}\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.475\n",
      "[ 0.01878417 -0.01260595  0.04084862 -0.03266317  0.01502334 -0.03667437\n",
      "  0.01255934  0.04815065 -0.02206419  0.02479605  0.06899284  0.0193719\n",
      " -0.01988549 -0.0087049   0.04605863  0.05793382  0.061218   -0.04720391\n",
      "  0.06070375 -0.01610907 -0.03484607]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "train=[]\n",
    "test=[]\n",
    "\n",
    "#数据读入\n",
    "with open('hw3_train.dat') as file:\n",
    "    for i in file.readlines():\n",
    "        train.append([1]+list(map(float,i.strip().split(' '))))\n",
    "\n",
    "with open('hw3_test.dat') as file:\n",
    "    for i in file.readlines():\n",
    "        test.append([1]+list(map(float,i.strip().split(' '))))\n",
    "        \n",
    "train=np.array(train)\n",
    "test=np.array(test)\n",
    "\n",
    "#定义函数\n",
    "def f(y,w,x):\n",
    "    temp=y*w.dot(x)\n",
    "    return (-y*x)/(np.exp(temp)+1)\n",
    "\n",
    "def sig(w,x):\n",
    "    return 1/(math.exp(-w.dot(x))+1)\n",
    "\n",
    "#数据维度\n",
    "m=train.shape[1]-1\n",
    "#数据组数\n",
    "n=train.shape[0]\n",
    "\n",
    "w=np.zeros(m)\n",
    "k=0.001\n",
    "\n",
    "for i in range(2000):\n",
    "    s=np.zeros(m)\n",
    "    for j in range(n):\n",
    "        s+=f(train[j][-1],w,train[j][:-1])\n",
    "    s=s/n\n",
    "    w-=k*s\n",
    "\n",
    "#计算Xw\n",
    "r1=test[:,:-1].dot(w)\n",
    "#计算sign(Xw)\n",
    "r2=np.sign(r1)\n",
    "#求出误差\n",
    "print((r2!=test[:,-1]).sum()/test.shape[0])\n",
    "print(w)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "所以$E_{out}$为0.475"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 19\n",
    "\n",
    "取$\\eta=0.01$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.22\n",
      "[-0.00385379 -0.18914564  0.26625908 -0.35356593  0.04088776 -0.3794296\n",
      "  0.01982783  0.33391527 -0.26386754  0.13489328  0.4914191   0.08726107\n",
      " -0.25537728 -0.16291797  0.30073678  0.40014954  0.43218808 -0.46227968\n",
      "  0.43230193 -0.20786372 -0.36936337]\n"
     ]
    }
   ],
   "source": [
    "w=np.zeros(m)\n",
    "k=0.01\n",
    "\n",
    "for i in range(2000):\n",
    "    s=np.zeros(m)\n",
    "    for j in range(n):\n",
    "        s+=f(train[j][-1],w,train[j][:-1])\n",
    "    s=s/n\n",
    "    w-=k*s\n",
    "\n",
    "#计算Xw\n",
    "r1=test[:,:-1].dot(w)\n",
    "#计算sign(Xw)\n",
    "r2=np.sign(r1)\n",
    "#求出误差\n",
    "print((r2!=test[:,-1]).sum()/test.shape[0])\n",
    "print(w)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "所以$E_{out}$为0.22"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 20\n",
    "\n",
    "使用随机梯度下降，只要对之前的式子稍作修改即可。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.477\n",
      "[ 0.03102133 -0.00824134  0.05216681 -0.02260935  0.02279948 -0.03096168\n",
      "  0.01339008  0.05210572 -0.00518606  0.02777817  0.07374788  0.03573774\n",
      " -0.01246505 -0.00478996  0.04042543  0.06805246  0.07474674 -0.04407857\n",
      "  0.06226404 -0.01511765 -0.03630915]\n"
     ]
    }
   ],
   "source": [
    "w=np.zeros(m)\n",
    "k=0.001\n",
    "\n",
    "for i in range(2000):\n",
    "    j=np.random.choice(n)\n",
    "    s=f(train[j][-1],w,train[j][:-1])\n",
    "    w-=k*s\n",
    "\n",
    "#计算Xw\n",
    "r1=test[:,:-1].dot(w)\n",
    "#计算sign(Xw)\n",
    "r2=np.sign(r1)\n",
    "#求出误差\n",
    "print((r2!=test[:,-1]).sum()/test.shape[0])\n",
    "print(w)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "所以$E_{out}$为0.477"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "以下两题为附加题\n",
    "\n",
    "#### Problem 21\n",
    "\n",
    "先回顾题目\n",
    "$$\n",
    "h = (h(x_1), h(x_2),... , h(x_N))\\\\\n",
    "y = (y_1,y_2, ... , y_N)\\\\\n",
    "RMSE(h) =\\sqrt{\\frac 1N \\sum_{i=1}^{N}(y_n − h(x_n))^2}\n",
    "$$\n",
    "题目问的是要计算$h^Ty$，至少需要调用几次$RMSE(h) $，注意这里只知道$h$。首先肯定要求$y$，因为有$N$个未知数，所以第一感觉是要调用$N$次，但是$N=1$时就不成立，因为有平方项。所以推测调用$N$次不行，接下来证明这个结论，注意这里$h,y$为行向量，所以和一般看到的形式略有不同。\n",
    "\n",
    "对$RMSE(h) $进行改写\n",
    "$$\n",
    "\\begin{aligned}\n",
    "RMSE(h) &=\\sqrt{\\frac 1N \\sum_{i=1}^{N}(y_n − h(x_n))^2}\n",
    "\\\\&=\\sqrt{\\frac 1N ||h-y||^2}\n",
    "\\\\&=\\sqrt{\\frac 1N (h-y)(h-y)^T}\n",
    "\\\\&=\\sqrt{\\frac 1N (yy^T-2hy^T+hh^T)}\n",
    "\\end{aligned}\n",
    "$$\n",
    "两边平方移项可得\n",
    "$$\n",
    "yy^T-2hy^T+hh^T=N\\times (RMSE(h))^2\n",
    "$$\n",
    "现在对两个不同的$h_1,h_2$调用$RMSE(h) $\n",
    "$$\n",
    "yy^T-2h_1y^T+h_1h_1^T=N\\times (RMSE(h_1))^2\\\\\n",
    "yy^T-2h_2y^T+h_2h_2^T=N\\times (RMSE(h_2))^2\\\\\n",
    "两式相减可得\\\\\n",
    "2(h_2-h_1)y^T=N\\times (RMSE(h_1))^2-N\\times (RMSE(h_2))^2\n",
    "$$\n",
    "这样就得到了一个线性方程。现在对$h_1,h_2,...,h_{k}$分别调用$RMSE(h) $，计算$RMSE(h_i)-RMSE(h_1)$，其中$(i=2,...k)$，根据之前所述可以得到$k-1$个线性方程组，有如下形式\n",
    "$$\n",
    "M_1y^T=M_2\\\\\n",
    "M_1\\in R^{(k-1)\\times N},y^T\\in R^N,M_2\\in R^{k-1}\n",
    "$$\n",
    "由线性代数知识我们知道\n",
    "$$\n",
    "当r(M_1)<N时，该方程有无穷多组解，当r(M_1)=N时，该方程有唯一解，其中r(M)表示矩阵M的秩。\n",
    "$$\n",
    "所以当$k=N$时，$r(M_1)\\le k-1<N$，此时有无穷多组解，所以无法确定$y = (y_1,y_2, ... , y_N)$；当$k=N+1$时，$r(M_1)$可能为$N$，若$r(M_1)=N$，此时可以唯一解出$y = (y_1,y_2, ... , y_N)$。\n",
    "\n",
    "综上，至少需要调用$N+1$次$RMSE(h) $。\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Problem 22\n",
    "\n",
    "为方便叙述，这里做以下记号，注意这里为列向量，和上题有所不同\n",
    "$$\n",
    "h= \\left(\n",
    " \\begin{matrix}\n",
    "   h_1(x_1) & h_2(x_1) & ... & h_K(x_1)\\\\\n",
    "   h_1(x_2) & h_2(x_2) & ... & h_K(x_2) \\\\\n",
    "     ...  &  ...  & ... &  ...  \\\\\n",
    " h_1(x_N) & h_2(x_N) & ... & h_K(x_N)\n",
    "  \\end{matrix}\n",
    "  \\right) \\in R^{N\\times K}\\\\\n",
    "  w=(w_1,w_2,...,w_K)^T\\in R^{ K}\\\\\n",
    "  y= (y_1,y_2, ... , y_N)^T\n",
    "$$\n",
    "那么$RMSE(H) $可以表示为\n",
    "$$\n",
    "RMSE(H) =\\sqrt{\\frac 1N ||y-hw||^2}=\\sqrt{\\frac 1N (y-hw)^T(y-hw)}\n",
    "$$\n",
    "由线性回归的推导我们知道最小化$RMSE(H)$的$w$满足以下条件\n",
    "$$\n",
    "X^TXw=X^Ty\n",
    "$$\n",
    "$X$已知，$y$未知，所以求出$y$即可，由上一题我们知道至少调用$N+1$次$RMSE(H)$就可以求出$y$，所以这题也至少需要调用$N+1$次。"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
